diff --git a/mcstas-comps/examples/ESS/ESS_MCPL_generate_filtered/ESS_MCPL_generate_filtered.instr b/mcstas-comps/examples/ESS/ESS_MCPL_generate_filtered/ESS_MCPL_generate_filtered.instr new file mode 100644 index 000000000..c66a83c9b --- /dev/null +++ b/mcstas-comps/examples/ESS/ESS_MCPL_generate_filtered/ESS_MCPL_generate_filtered.instr @@ -0,0 +1,241 @@ +/******************************************************************************* +* Instrument: ESS_MCPL_generate_filtered +* +* %I +* Written by: Peter Willendrup +* Date: 2016-09-26 +* Origin: ESS +* %INSTRUMENT_SITE: ESS +* +* Instrument meant to 'filter' legacy target division MCPL files. +* +* %D +* Instrument meant to 'filter' legacy target division MCPL files (i.e. S2.mcpl.gz -> S2_filtered.mcpl.gz) for use with +* e.g. the ESS_butterfly_MCPL_test instrument and filter=1. +* +* To function, the instrument needs a local copy of the MCPL description you want to work with, e.g. from +* https://public.esss.dk/users/willend/MCPL/ +* +* The instrument assumes that binary MCPL datasets are available in . named [sector][beamline].mcpl.gz, i.e. W8.mcpl.gz. +* +* %P +* sector: [str] Defines the 'sector' of your instrument position. Valid values are "N","S","E" and "W" +* beamline: [1] Defines the 'beamline number' of your instrument position. Valid values are 1..10 or 1..11 depending on sector +* Lmin: [AA] Minimum wavelength simulated +* Lmax: [AA] Maximum wavelength simulated +* c_performance: [1] Cold brilliance scalar performance multiplicator c_performance > 0 +* t_performance: [1] Thermal brilliance scalar performance multiplicator t_performance > 0 +* index: [1] Target index for source focusing. Defaults to illuminate the "cold collimated" brilliance monitor, thereby suppressing "dist" +* dist: [m] Distance from origin to focusing rectangle; at (0,0,dist) - alternatively use target_index +* cold: [1] Defines the statistical fraction of events emitted from the cold part of the moderator +* Yheight: [m] Defines the moderator height. Valid values are 0.03 m and 0.06 m +* delta: [m] Parameter that allows to scan "collimator" position +* allmons: [1] Flag to include also debugging-monitor output +* +* %L +* +* Benchmarking website available at http://ess_butterfly.mcstas.org +* %E +*******************************************************************************/ +DEFINE INSTRUMENT ESS_MCPL_generate_filtered(string sector="S", + int beamline=2,Lmin=0.2,Lmax=20,c_performance=1,t_performance=1, + int index=0,dist=2,cold=0.5,Yheight=0.03,delta=0,int allmons=0) + + +DECLARE %{ + double calcAlpha(double length, double radius) { + // calculate angle of arm after curved guide + return RAD2DEG * length/radius; + } + + double calcX(double length, double radius) { + // calculate position and angle of arm after curved guide + double alpha = DEG2RAD * calcAlpha(length, radius); + return radius*(1.0-cos(alpha)); + } + + double calcZ(double length, double radius) { + // calculate position and angle of arm after curved guide + double alpha = DEG2RAD * calcAlpha(length, radius); + return radius*sin(alpha); + } + + double XW, YH; + char options1[256],options2[256],options3[256],options4[256]; + char srcdef[128]; + double WidthC=0.072,WidthT=0.108; + double lambdamin, lambdamax; + double TCollmin; + double TCollmax; + #pragma acc declare create(TCollmin,TCollmax) + double EminTh=20, EmaxTh=100, EminC=0, EmaxC=20; + double weightmultiplier; + #pragma acc declare create(EminTh,EmaxTh,EminC,EmaxC,weightmultiplier) + /* 10 beamlines in sector N and E - plus one location added for drawing */ + double iBeamlinesN[] = { 30.0, 36.0, 42.0, 48.0, 54.0, 60.0, 66.0, 72.0, 78.0, 84.0, 90.0}; + double iBeamlinesE[] = {-30.0, -36.0, -42.0, -48.0, -54.0, -60.0, -66.0, -72.0, -78.0, -84.0, -90.0}; + /* 11 beamlines in sector S and W - plus one location added for drawing */ + double iBeamlinesW[] = { 150.0, 144.7, 138.0, 132.7, 126.0, 120.7, 114.0, 108.7, 102.0, 96.7, 90.0, 84.0}; + double iBeamlinesS[] = {-150.0, -144.7, -138.0, -132.7, -126.0, -120.7, -114.0, -108.7, -102.0, -96.7, -90.0, -84.0}; + double* iBeamlines; + double ANGLE; + double DeltaX,DeltaZ; + char MCPLfile[128]; + char FILTERfile[128]; + int Divisions; +%} + +USERVARS %{ + int IsCold; + double SrcX; + double SrcY; + double SrcZ; + double E_min; + double E_max; + double Eneutron; + double T0; + double L0; + /* particle data in original MCNP / TCS units */ + double xorig; + double yorig; + double zorig; + double vxorig; + double vyorig; + double vzorig; + double torig; +%} + +INITIALIZE +%{ + lambdamin=Lmin; + lambdamax=Lmax; + XW=1.05*(WidthC+2*WidthT); + YH=1.05*Yheight; + sprintf(options1,"user1 bins=201 limits=[-%g,%g]",XW/2,XW/2); + sprintf(options4,"user1 bins=201 limits=[-%g,%g]",YH/2,YH/2); + sprintf(options2,"user1 bins=201 limits=[-%g,%g], user2 bins=201 limits=[-%g,%g]",XW/2,XW/2,YH/2,YH/2); + sprintf(options3,"user1 bins=201 limits=[-%g,%g], user2 bins=201 limits=[-%g,%g]",1.05*(WidthC/2),1.05*(WidthC/2),1.05*Yheight/2,1.05*Yheight/2); + sprintf(srcdef,"2015"); + if (beamline==1) { + TCollmin=0; + TCollmax=0.058; + } else if (beamline==2) { + TCollmin=0; + TCollmax=0.06; + } + else { + TCollmin=0.011; + TCollmax=0.071; + } + #pragma acc update device(TCollmin,TCollmax) + if (strcasestr(sector,"N")) { + iBeamlines=iBeamlinesN; + DeltaX=-0.0585; DeltaZ=0.0925; + } else if (strcasestr(sector,"W")) { + iBeamlines=iBeamlinesW; + DeltaX=0.0585; DeltaZ=0.0925; + } else if (strcasestr(sector,"S")) { + iBeamlines=iBeamlinesS; + DeltaX=0.0585; DeltaZ=-0.0925; + } else if (strcasestr(sector,"E")) { + iBeamlines=iBeamlinesE; + DeltaX=-0.0585; DeltaZ=-0.0925; + } + ANGLE=iBeamlines[beamline-1]-90; + sprintf(MCPLfile,"%s%i.mcpl.gz",sector,beamline); + sprintf(FILTERfile,"%s%i_filtered.mcpl.gz",sector,beamline); + printf("MCPLfile is %s\n",MCPLfile); + printf("Filtered MCPLfile is %s\n",FILTERfile); + + weightmultiplier=1.56e16/1e5; /* (ESS protons/s) / MCNP nps */ + #pragma acc update device(EminTh,EmaxTh,EminC,EmaxC,weightmultiplier) +%} + +TRACE + +COMPONENT Origin = Progress_bar() +AT (0, 0, 0) ABSOLUTE + +/* read neutrons from an mcpl file*/ + +COMPONENT vinROT2 = Arm() +AT(0,0,0) RELATIVE PREVIOUS + ROTATED (0,-90,0) RELATIVE PREVIOUS + +COMPONENT vinROT1 = Arm() +AT(0,0,0) RELATIVE PREVIOUS + ROTATED (-90,0,0) RELATIVE PREVIOUS + + +COMPONENT vin = MCPL_input_once(filename=MCPLfile,verbose=1) +AT(0,0,-0.137) RELATIVE PREVIOUS +EXTEND %{ + SCATTER; + xorig=x; yorig=y; zorig=z; + vxorig=vx; vyorig=vy; vzorig=vz; + torig=t; +%} + +COMPONENT Sphere1 = PSD_monitor_4PI(filename="nonrotated", radius=2.2,restore_neutron=1, nowritefile=!allmons) +WHEN (allmons) AT (0,0,0) RELATIVE PREVIOUS + +COMPONENT Source = ESS_butterfly(sector=sector,beamline=beamline,Lmin=Lmin,Lmax=Lmax,c_performance=c_performance,t_performance=t_performance,dist=dist,target_index=index,cold_frac=cold, yheight=Yheight, + focus_xw=0.12, focus_yh=0.12) + WHEN (0==1) AT (DeltaX,0,DeltaZ) ABSOLUTE + ROTATED (0, ANGLE, 0) ABSOLUTE + +COMPONENT Sphere0 = PSD_monitor_4PI(filename="rotated", radius=2.2,restore_neutron=1, nowritefile=!allmons) +WHEN (allmons) AT (0,0,0) RELATIVE Source + +COMPONENT Focus_cut=Shape(xwidth=0.01,yheight=0.01) +AT(0,0,2) RELATIVE Source +EXTEND %{ + ALLOW_BACKPROP; + PROP_Z0; + SCATTER; + if (fabs(x)>0.06 || fabs(y)>0.06) { + ABSORB; + } +%} + +COMPONENT ATSource=PSD_monitor(xwidth=0.31,yheight=0.31, nowritefile=!allmons) + WHEN (allmons) AT (0,0,0.08) RELATIVE Source + +COMPONENT BackTrace = Shape(xwidth=0.3,yheight=0.3) +AT (0,0,0.08) RELATIVE Source +EXTEND %{ + /* Propagate back to a small rectangle in front of moderators */ + SCATTER; + double myL = (2*PI/V2K)/sqrt(vx*vx + vy*vy + vz*vz); + ALLOW_BACKPROP; + PROP_Z0; + SCATTER; + /* Remove neutrons that are not from around the moderators */ + if (fabs(x)>0.12 || fabs(y)>0.03) { + ABSORB; + } + if ( myL < INSTRUMENT_GETPAR(Lmin) || myL > INSTRUMENT_GETPAR(Lmax) ) { + ABSORB; + } +%} + +COMPONENT POSTATSource=PSD_monitor(xwidth=0.31,yheight=0.31) +AT (0,0,0) RELATIVE BackTrace + +COMPONENT Restore = Arm() + AT (0,0,0) RELATIVE vin +EXTEND %{ + x=xorig; y=yorig; z=zorig; + vx=vxorig; vy=vyorig; vz=vzorig; + t=torig; +%} + +COMPONENT vout = MCPL_output(filename=FILTERfile) + AT (0,0,0) RELATIVE vin + + +FINALLY +%{ +%} + +END diff --git a/mcstas-comps/examples/ESS/ESS_MCPL_generate_filtered/README.md b/mcstas-comps/examples/ESS/ESS_MCPL_generate_filtered/README.md new file mode 100644 index 000000000..bb2244438 --- /dev/null +++ b/mcstas-comps/examples/ESS/ESS_MCPL_generate_filtered/README.md @@ -0,0 +1,49 @@ +# The `ESS_MCPL_generate_filtered` Instrument + +*McStas: Instrument meant to 'filter' legacy target division MCPL files.* + +## Identification + +- **Site:** ESS +- **Author:** Peter Willendrup +- **Origin:** ESS +- **Date:** 2016-09-26 + +## Description + +```text +Instrument meant to 'filter' legacy target division MCPL files (i.e. S2.mcpl.gz -> S2_filtered.mcpl.gz) for use with +e.g. the ESS_butterfly_MCPL_test instrument and filter=1. + +To function, the instrument needs a local copy of the MCPL description you want to work with, e.g. from +https://public.esss.dk/users/willend/MCPL/ + +The instrument assumes that binary MCPL datasets are available in . named [sector][beamline].mcpl.gz, i.e. W8.mcpl.gz. +``` + +## Input parameters + +Parameters in **boldface** are required; the others are optional. + +| Name | Unit | Description | Default | +|------|------|-------------|---------| +| sector | str | Defines the 'sector' of your instrument position. Valid values are "N","S","E" and "W" | "S" | +| beamline | 1 | Defines the 'beamline number' of your instrument position. Valid values are 1..10 or 1..11 depending on sector | 2 | +| Lmin | AA | Minimum wavelength simulated | 0.2 | +| Lmax | AA | Maximum wavelength simulated | 20 | +| c_performance | 1 | Cold brilliance scalar performance multiplicator c_performance > 0 | 1 | +| t_performance | 1 | Thermal brilliance scalar performance multiplicator t_performance > 0 | 1 | +| index | 1 | Target index for source focusing. Defaults to illuminate the "cold collimated" brilliance monitor, thereby suppressing "dist" | 0 | +| dist | m | Distance from origin to focusing rectangle; at (0,0,dist) - alternatively use target_index | 2 | +| cold | 1 | Defines the statistical fraction of events emitted from the cold part of the moderator | 0.5 | +| Yheight | m | Defines the moderator height. Valid values are 0.03 m and 0.06 m | 0.03 | +| delta | m | Parameter that allows to scan "collimator" position | 0 | +| allmons | 1 | Flag to include also debugging-monitor output | 0 | + +## Links + +- [Source code](ESS_MCPL_generate_filtered.instr) for `ESS_MCPL_generate_filtered.instr`. +- +- Benchmarking website available at http://ess_butterfly.mcstas.org + +--- diff --git a/mcstas-comps/examples/ESS/ESS_butterfly_Guide_curved_test/ESS_butterfly_Guide_curved_test.instr b/mcstas-comps/examples/ESS/ESS_butterfly_Guide_curved_test/ESS_butterfly_Guide_curved_test.instr index 9b53ecb8e..14660aaf6 100644 --- a/mcstas-comps/examples/ESS/ESS_butterfly_Guide_curved_test/ESS_butterfly_Guide_curved_test.instr +++ b/mcstas-comps/examples/ESS/ESS_butterfly_Guide_curved_test/ESS_butterfly_Guide_curved_test.instr @@ -7,13 +7,17 @@ * Origin: ESS * %INSTRUMENT_SITE: ESS * -* Test instrument for the updated BF1 butterfly moderator design +* Test instrument for the updated BF1 butterfly moderator design including a simple guide. * * %D -* Test instrument for the updated BF1 butterfly moderator design. +* Test instrument for the updated BF1 butterfly moderator design including a simple guide with 10 m straight flight-path followed by a 50 m curved guide. * -* The below example gives a 50-50 (statistics-wise) cold/thermal beam at beamline N10. -* %Example: sector=N beamline=10 cold=0.5 Detector: Monitor2_xy1_I=4.4e+11 +* The below example gives a 50-50 (statistics-wise) cold/thermal beam at beamline S2. +* Example: sector=S beamline=2 cold=0.5 Detector: Monitor2_xy1_I=1.59e+11 +* +* Use the ESS_butterfly_MCPL_test instrument for direct comparison with MCPL-based source descriptions. +* +* %Example: sector=S beamline=2 cold=0.5 Detector: Monitor2_xy1_I=1.59e+11 * * * %P @@ -28,22 +32,17 @@ * cold: [1] Defines the statistical fraction of events emitted from the cold part of the moderator * Yheight: [m] Defines the moderator height. Valid values are 0.03 m and 0.06 m * delta: [m] Parameter that allows to scan "collimator" position -* thres: [ ] Weight-threshold, neutrons with higher weight are absorbed -* repeat: [ ] Number of MCPL file repetitions -* v_smear: [1] When repeating events, make a Gaussian MC choice within v_smear*V around particle velocity V -* pos_smear: [m] When repeating events, make a flat MC choice of position within pos_smear around particle starting position -* dir_smear: [deg] When repeating events, make a Gaussian MC choice of direction within dir_smear around particle direction * n_pulses: [1] Number of pulses to simulate +* allmons: [1] Flag to include also debugging-monitor output * * %L * * Benchmarking website available at http://ess_butterfly.mcstas.org * %E *******************************************************************************/ -DEFINE INSTRUMENT ESS_butterfly_Guide_curved_test(string sector="N", - int beamline=1,Lmin=0.2,Lmax=20,c_performance=1,t_performance=1, - int index=25,dist=0,cold=1,Yheight=0.03,delta=0,thres=0.003,repeat=1, - v_smear=0.1,pos_smear=0.01,dir_smear=0.01,int n_pulses=1) +DEFINE INSTRUMENT ESS_butterfly_Guide_curved_test(string sector="S", + int beamline=2,Lmin=0.2,Lmax=20,c_performance=1,t_performance=1, + int index=0,dist=2,cold=0.5,Yheight=0.03,delta=0,int n_pulses=1,int allmons=0) DECLARE %{ @@ -90,8 +89,8 @@ USERVARS %{ double SrcX; double SrcY; double SrcZ; - double Emin; - double Emax; + double E_min; + double E_max; double Eneutron; double T0; double L0; @@ -142,7 +141,6 @@ TRACE COMPONENT Origin = Progress_bar() AT (0, 0, 0) ABSOLUTE - COMPONENT vinROT2 = Arm() AT(0,0,0) RELATIVE PREVIOUS ROTATED (0,-90,0) RELATIVE PREVIOUS @@ -154,14 +152,13 @@ AT(0,0,0) RELATIVE PREVIOUS COMPONENT vin = Arm() AT(0,0,0) RELATIVE PREVIOUS -/* Focusing for this use of the source is a little unphysical: 1x1cm @ 1m ~ 1e-4 steradian. To be useful in a "proper" instrument, you should of course illuminate your beamport fully!*/ COMPONENT Source = ESS_butterfly(sector=sector,beamline=beamline,Lmin=Lmin,Lmax=Lmax,c_performance=c_performance,t_performance=t_performance,dist=dist,target_index=index,cold_frac=cold, yheight=Yheight, focus_xw=0.12, focus_yh=0.12, n_pulses=n_pulses) AT (DeltaX,0,DeltaZ) ABSOLUTE ROTATED (0, ANGLE, 0) ABSOLUTE -COMPONENT Sphere0 = PSD_monitor_4PI(filename="rotated", radius=2.2,restore_neutron=1) -AT (0,0,0) RELATIVE Source +COMPONENT Sphere0 = PSD_monitor_4PI(filename="rotated", radius=2.2,restore_neutron=1, nowritefile=!allmons) +WHEN (allmons) AT (0,0,0) RELATIVE Source COMPONENT BackTrace = Shape(xwidth=0.3,yheight=0.3) AT (0,0,0.08) RELATIVE Source @@ -172,16 +169,18 @@ EXTEND %{ SrcX=x;SrcY=y;SrcZ=z; Eneutron=VS2E*(vx*vx + vy*vy + vz*vz); if (Eneutron>EminTh) { - Emin=EminC;Emax=EmaxC; + E_min=EminC;E_max=EmaxC; IsCold=0; } else { - Emin=EminTh;Emax=EmaxTh; + E_min=EminTh;E_max=EmaxTh; IsCold=1; } T0=t; L0=myL; %} +COMPONENT POSTATSource=PSD_monitor(xwidth=0.31,yheight=0.31) +AT (0,0,0) RELATIVE BackTrace /* These arms are just to ensure we get a good view of the monolith */ COMPONENT Arm1 = Arm() @@ -190,49 +189,48 @@ COMPONENT Arm1 = Arm() COMPONENT Arm2 = Arm() AT (0,0,3.5) RELATIVE Source -COMPONENT AutoTOFL0 = Monitor_nD(xwidth=XW, yheight=YH, user1="T0", user2="L0", options="user1 limits=[0 5e-3] bins=51, user2 limits=[0.1 20] bins=41", restore_neutron=1) +COMPONENT AutoTOFL0 = Monitor_nD(xwidth=XW, yheight=YH, user1="T0", user2="L0", options="user1 limits=[0 5e-3] bins=51, user2 limits=[0.1 20] bins=41", restore_neutron=1,filename="AutoTOFL0") AT (0, 0, 0.001) RELATIVE BackTrace -COMPONENT AutoTOF0 = Monitor_nD(xwidth=XW, yheight=YH, user1="T0", options="user1 limits=[0 5e-3] bins=51", restore_neutron=1) +COMPONENT AutoTOF0 = Monitor_nD(xwidth=XW, yheight=YH, user1="T0", options="user1 limits=[0 5e-3] bins=51", restore_neutron=1,filename="AutoTOF0") AT (0, 0, 0.001) RELATIVE PREVIOUS -COMPONENT AutoL0 = Monitor_nD(xwidth=XW, yheight=YH, user1="L0", options="user1 limits=[0.1 20] bins=41", restore_neutron=1) +COMPONENT AutoL0 = Monitor_nD(xwidth=XW, yheight=YH, user1="L0", options="user1 limits=[0.1 20] bins=41", restore_neutron=1,filename="AutoL0") AT (0, 0, 0.001) RELATIVE PREVIOUS -COMPONENT PSD0= Monitor_nD(filename="flat",xwidth=0.2,yheight=0.2,user1="SrcX",user2="SrcY",options="user1 limits=[-0.1 0.1] bins=90, user2 limits=[-0.1 0.1] bins=90,", restore_neutron=1) +COMPONENT PSD0= Monitor_nD(filename="flat",xwidth=0.2,yheight=0.2,user1="SrcX",user2="SrcY",options="user1 limits=[-0.1 0.1] bins=90, user2 limits=[-0.1 0.1] bins=90,", restore_neutron=1,filename="PSD0") AT (0,0,0.001) RELATIVE PREVIOUS -COMPONENT PSD1=Monitor_nD(filename="flatC",xwidth=0.2,yheight=0.2,user1="SrcX",user2="SrcY",options="user1 limits=[-0.1 0.1] bins=90, user2 limits=[-0.1 0.1] bins=90,", restore_neutron=1) +COMPONENT PSD1=Monitor_nD(filename="flatC",xwidth=0.2,yheight=0.2,user1="SrcX",user2="SrcY",options="user1 limits=[-0.1 0.1] bins=90, user2 limits=[-0.1 0.1] bins=90,", restore_neutron=1,filename="PSD1") WHEN (Eneutron=EminTh) AT (0,0,0.001) RELATIVE PREVIOUS -COMPONENT MonND1 = Monitor_nD(xwidth=XW, yheight=YH, user1="SrcX", username1="Horizontal position / [m]", options=options1, restore_neutron=1) +COMPONENT MonND1 = Monitor_nD(xwidth=XW, yheight=YH, user1="SrcX", username1="Horizontal position / [m]", options=options1, restore_neutron=1,filename="MonND1") AT (0, 0, 1) RELATIVE Source -COMPONENT CWidth = Monitor_nD(xwidth=XW, yheight=YH, user1="SrcX", username1="Horizontal position / [m]", options=options1, restore_neutron=1) +COMPONENT CWidth = Monitor_nD(xwidth=XW, yheight=YH, user1="SrcX", username1="Horizontal position / [m]", options=options1, restore_neutron=1,filename="CWidth") WHEN(Eneutron<=EmaxC && Eneutron>=EminC) AT (0, 0, 1) RELATIVE Source -COMPONENT TWidth = Monitor_nD(xwidth=XW, yheight=YH, user1="SrcX", username1="Horizontal position / [m]", options=options1, restore_neutron=1) +COMPONENT TWidth = Monitor_nD(xwidth=XW, yheight=YH, user1="SrcX", username1="Horizontal position / [m]", options=options1, restore_neutron=1,filename="TWidth") WHEN(Eneutron<=EmaxTh && Eneutron>=EminTh) AT (0, 0, 1) RELATIVE Source -COMPONENT MonND2 = Monitor_nD(xwidth=XW, yheight=YH, user1="SrcY", username1="Vertical position COLD / [m]", options=options4, restore_neutron=1) +COMPONENT MonND2 = Monitor_nD(xwidth=XW, yheight=YH, user1="SrcY", username1="Vertical position COLD / [m]", options=options4, restore_neutron=1,filename="MonND2") WHEN(IsCold) AT (0, 0, 1) RELATIVE Source -COMPONENT MonND2_2 = Monitor_nD(xwidth=XW, yheight=YH, user1="SrcY", username1="Vertical position THERMAL/ [m]", options=options4, restore_neutron=1) +COMPONENT MonND2_2 = Monitor_nD(xwidth=XW, yheight=YH, user1="SrcY", username1="Vertical position THERMAL/ [m]", options=options4, restore_neutron=1,filename="MonND2_2") WHEN(!IsCold) AT (0, 0, 1) RELATIVE Source -COMPONENT MonND3 = Monitor_nD(xwidth=XW, yheight=YH, user1="SrcX", username1="Horizontal position / [m]", user2="SrcY",username2="Vertical position / [m]", options=options2, restore_neutron=1) +COMPONENT MonND3 = Monitor_nD(xwidth=XW, yheight=YH, user1="SrcX", username1="Horizontal position / [m]", user2="SrcY",username2="Vertical position / [m]", options=options2, restore_neutron=1,filename="MonND3") AT (0, 0, 1) RELATIVE Source -COMPONENT MonND4 = Monitor_nD(xwidth=XW, yheight=YH, user1="SrcX", username1="Emission position / [m]", user2="SrcZ", username2="Z-component of position / [m]", options="user1 bins=201 limits=[-0.3,0.3], user2 bins=201 limits=[-0.3,0.3]", restore_neutron=1) +COMPONENT MonND4 = Monitor_nD(xwidth=XW, yheight=YH, user1="SrcX", username1="Emission position / [m]", user2="SrcZ", username2="Z-component of position / [m]", options="user1 bins=201 limits=[-0.3,0.3], user2 bins=201 limits=[-0.3,0.3]", restore_neutron=1,filename="MonND4") AT (0, 0, 1) RELATIVE Source - -COMPONENT AutoTOFL = Monitor_nD(xwidth=0.1, yheight=0.1, options="tof limits=[0 15e-3] bins=51, lambda limits=[0.1 20] bins=41", restore_neutron=1) +COMPONENT AutoTOFL = Monitor_nD(xwidth=0.1, yheight=0.1, options="tof limits=[0 15e-3] bins=51, lambda limits=[0.1 20] bins=41", restore_neutron=1,filename="AutoTOFL") AT (0, 0, 2) RELATIVE Source @@ -241,8 +239,8 @@ COMPONENT AutoTOFL = Monitor_nD(xwidth=0.1, yheight=0.1, options="tof limits=[0 COMPONENT BrillmonCOLD = Brilliance_monitor( nlam = 101, nt = 101, filename = "brillCOLD", t_0 = -1000, t_1 =4e4, lambda_0 = lambdamin, lambda_1 = lambdamax, - Freq =14, toflambda=1 ,tofcuts=0, srcarea=(100*0.072*100*Yheight), restore_neutron=1,source_dist=2,xwidth=0.1,yheight=0.1) -WHEN(IsCold) AT (0, 0, 2) RELATIVE Source + Freq =14, toflambda=1 ,tofcuts=0, srcarea=(100*0.072*100*Yheight), restore_neutron=1,source_dist=2,xwidth=0.1,yheight=0.1,nowritefile=!allmons) +WHEN(allmons && IsCold) AT (0, 0, 2) RELATIVE Source /* Measures "collimated" brilliance of the cold source over fixed 6 cm wide area x central part vertically. */ /* Used for calibration of performance wrt. MCNP BF1 output, see http://ess_butterfly.mcstas.org */ @@ -256,8 +254,8 @@ COMPONENT BrillmonCOLD_COLL = Brilliance_monitor( COMPONENT BrillmonTHRM = Brilliance_monitor( nlam = 101, nt = 101, filename = "brillTHRM", t_0 = -1000, t_1 =4e4, lambda_0 = lambdamin, lambda_1 = lambdamax, - Freq =14, toflambda=1,tofcuts=0, srcarea=(100*0.108*100*Yheight), restore_neutron=1,source_dist=2,xwidth=0.1,yheight=0.1) - WHEN (!IsCold) AT (0, 0, 2) RELATIVE Source + Freq =14, toflambda=1,tofcuts=0, srcarea=(100*0.108*100*Yheight), restore_neutron=1,source_dist=2,xwidth=0.1,yheight=0.1,nowritefile=!allmons) + WHEN (allmons && !IsCold) AT (0, 0, 2) RELATIVE Source /* Measures "collimated" brilliance of the thermal source over fixed 6 cm wide area (or smaller at beamlines no. 1,2) x central part vertically. */ /* Used for calibration of performance wrt. MCNP BF1 output, see http://ess_butterfly.mcstas.org */ @@ -273,12 +271,20 @@ COMPONENT PSD0x= Monitor_nD(filename="flat_x",xwidth=0.1,yheight=0.1,user1="SrcX COMPONENT PSD1x=Monitor_nD(filename="flatC_x",xwidth=0.1,yheight=0.1,user1="SrcX",user2="SrcY",options="user1 limits=[-0.08 0.08] bins=90, user2 limits=[-0.08 0.08] bins=90,", restore_neutron=1) WHEN (Eneutron=EminTh) AT (0,0,0.001) RELATIVE PREVIOUS +COMPONENT Flux_incoming=Monitor_nD(xwidth=0.05,yheight=0.1,filename="Flux_incoming",options="x, y pr cm2", bins=100) + AT (0, 0, 2) RELATIVE Source + +COMPONENT GuideStraight=Guide(w1 = 0.05, h1 = 0.1, l = 10) + AT (0,0,0) RELATIVE PREVIOUS + +COMPONENT Flux_mid=Monitor_nD(xwidth=0.05,yheight=0.1,filename="Flux_mid",options="x, y pr cm2", bins=100) + AT (0, 0, 10) RELATIVE GuideStraight COMPONENT GuideR = Guide_curved( - w1 = 0.1, h1 = 0.1, l = 50, curvature=3000) + w1 = 0.05, h1 = 0.1, l = 50, curvature=3000) AT (0,0,0.001) RELATIVE PREVIOUS COMPONENT RArm=Arm() @@ -290,9 +296,32 @@ COMPONENT GuideR = Guide_curved( ROTATED (0, calcAlpha(50,3000),0) RELATIVE GuideR COMPONENT Monitor2_xy1 = Monitor_nD( - options = "x limits=[-0.06 0.06] bins=51, y limits=[-0.06 0.06] bins=51,", xwidth = 0.12, yheight = 0.12) + options = "x limits=[-0.06 0.06] bins=51, y limits=[-0.06 0.06] bins=51,", xwidth = 0.12, yheight = 0.12,filename="Monitor2_xy1") AT (0, 0, 0.05) RELATIVE RArm +COMPONENT Monitor_X = Monitor_nD( + options = "x limits=[-0.06 0.06] bins=51", xwidth = 0.12, yheight = 0.12,filename="Monitor_X") + AT (0, 0, 0) RELATIVE Monitor2_xy1 + +COMPONENT Monitor_Y = Monitor_nD( + options = "y limits=[-0.06 0.06] bins=51", xwidth = 0.12, yheight = 0.12,filename="Monitor_Y") + AT (0, 0, 0) RELATIVE Monitor2_xy1 + +COMPONENT Monitor_divH = Monitor_nD( + options = "hdiv limits=[-1.5 1.5] bins=51", xwidth = 0.12, yheight = 0.12,filename="Monitor_divH") + AT (0, 0, 0) RELATIVE Monitor2_xy1 + +COMPONENT Monitor_divV = Monitor_nD( + options = "vdiv limits=[-1.5 1.5] bins=51", xwidth = 0.12, yheight = 0.12,filename="Monitor_divV") + AT (0, 0, 0) RELATIVE Monitor2_xy1 + +COMPONENT Monitor_t = Monitor_nD( + options = "t limits=[0 350e-3] bins=351", xwidth = 0.12, yheight = 0.12,filename="Monitor_t") + AT (0, 0, 0) RELATIVE Monitor2_xy1 + +COMPONENT Flux_end=Monitor_nD(xwidth=0.05,yheight=0.1,filename="Flux_end",options="x, y pr cm2", bins=100) + AT (0, 0, 0) RELATIVE Monitor2_xy1 + /* /\* Uncomment these helper-arms to view "full" monolith *\/ */ COMPONENT DummyArm1 = Arm() diff --git a/mcstas-comps/examples/ESS/ESS_butterfly_Guide_curved_test/README.md b/mcstas-comps/examples/ESS/ESS_butterfly_Guide_curved_test/README.md index 397bf0565..35de80218 100644 --- a/mcstas-comps/examples/ESS/ESS_butterfly_Guide_curved_test/README.md +++ b/mcstas-comps/examples/ESS/ESS_butterfly_Guide_curved_test/README.md @@ -1,6 +1,6 @@ # The `ESS_butterfly_Guide_curved_test` Instrument -*McStas: Test instrument for the updated BF1 butterfly moderator design* +*McStas: Test instrument for the updated BF1 butterfly moderator design including a simple guide.* ## Identification @@ -12,9 +12,12 @@ ## Description ```text -Test instrument for the updated BF1 butterfly moderator design. +Test instrument for the updated BF1 butterfly moderator design including a simple guide with 10 m straight flight-path followed by a 50 m curved guide. -The below example gives a 50-50 (statistics-wise) cold/thermal beam at beamline N10. +The below example gives a 50-50 (statistics-wise) cold/thermal beam at beamline S2. +Example: sector=S beamline=2 cold=0.5 Detector: Monitor2_xy1_I=1.59e+11 + +Use the ESS_butterfly_MCPL_test instrument for direct comparison with MCPL-based source descriptions. ``` ## Input parameters @@ -23,23 +26,19 @@ Parameters in **boldface** are required; the others are optional. | Name | Unit | Description | Default | |------|------|-------------|---------| -| sector | str | Defines the 'sector' of your instrument position. Valid values are "N","S","E" and "W" | "N" | -| beamline | 1 | Defines the 'beamline number' of your instrument position. Valid values are 1..10 or 1..11 depending on sector | 1 | +| sector | str | Defines the 'sector' of your instrument position. Valid values are "N","S","E" and "W" | "S" | +| beamline | 1 | Defines the 'beamline number' of your instrument position. Valid values are 1..10 or 1..11 depending on sector | 2 | | Lmin | AA | Minimum wavelength simulated | 0.2 | | Lmax | AA | Maximum wavelength simulated | 20 | | c_performance | 1 | Cold brilliance scalar performance multiplicator c_performance > 0 | 1 | | t_performance | 1 | Thermal brilliance scalar performance multiplicator t_performance > 0 | 1 | -| index | 1 | Target index for source focusing. Defaults to illuminate the "cold collimated" brilliance monitor, thereby suppressing "dist" | 25 | -| dist | m | Distance from origin to focusing rectangle; at (0,0,dist) - alternatively use target_index | 0 | -| cold | 1 | Defines the statistical fraction of events emitted from the cold part of the moderator | 1 | +| index | 1 | Target index for source focusing. Defaults to illuminate the "cold collimated" brilliance monitor, thereby suppressing "dist" | 0 | +| dist | m | Distance from origin to focusing rectangle; at (0,0,dist) - alternatively use target_index | 2 | +| cold | 1 | Defines the statistical fraction of events emitted from the cold part of the moderator | 0.5 | | Yheight | m | Defines the moderator height. Valid values are 0.03 m and 0.06 m | 0.03 | | delta | m | Parameter that allows to scan "collimator" position | 0 | -| thres | | Weight-threshold, neutrons with higher weight are absorbed | 0.003 | -| repeat | | Number of MCPL file repetitions | 1 | -| v_smear | 1 | When repeating events, make a Gaussian MC choice within v_smear*V around particle velocity V | 0.1 | -| pos_smear | m | When repeating events, make a flat MC choice of position within pos_smear around particle starting position | 0.01 | -| dir_smear | deg | When repeating events, make a Gaussian MC choice of direction within dir_smear around particle direction | 0.01 | | n_pulses | 1 | Number of pulses to simulate | 1 | +| allmons | 1 | Flag to include also debugging-monitor output | 0 | ## Links diff --git a/mcstas-comps/examples/ESS/ESS_butterfly_MCPL_test/ESS_butterfly_MCPL_test.instr b/mcstas-comps/examples/ESS/ESS_butterfly_MCPL_test/ESS_butterfly_MCPL_test.instr index 53fce4208..dd1a1e572 100644 --- a/mcstas-comps/examples/ESS/ESS_butterfly_MCPL_test/ESS_butterfly_MCPL_test.instr +++ b/mcstas-comps/examples/ESS/ESS_butterfly_MCPL_test/ESS_butterfly_MCPL_test.instr @@ -7,13 +7,19 @@ * Origin: ESS * %INSTRUMENT_SITE: ESS * -* Test instrument for the updated BF1 butterfly moderator design using MCPL input +* Test instrument for the ESS source using direct MCPL input. Includes a simple guide. * * %D -* Test instrument for the updated BF1 butterfly moderator design using MCPL input files. +* Test instrument for the updated BF1 butterfly moderator design using MCPL input files. The instrument includes a simple guide with 10 m straight flight-path followed by a 50 m curved guide. * -* The below example gives a 50-50 (statistics-wise) cold/thermal beam at beamline N10. -* Example: sector=N beamline=10 cold=0.5 +* The below example gives a 50-50 (statistics-wise) cold/thermal beam at beamline S2 using the base S2 MCPL file from +* https://public.esss.dk/users/willend/MCPL/S2.mcpl.gz +* Example: sector=S beamline=2 Detector: Monitor2_xy1_I=1.59e+11 +* and a 'preconditioned' version of this the file from https://public.esss.dk/users/willend/MCPL/S2_filtered.mcpl.gz +* (This file is fetched via a curl command at code-generation time) +* Example: sector=S beamline=2 filter=1 Detector: Monitor2_xy1_I=1.59e+11 +* %Example: sector=S beamline=2 filter=1 thres=0 Detector: Monitor2_xy1_I=1.59e+11 +* %Example: sector=S beamline=2 filter=1 thres=4e8 Detector: Monitor2_xy1_I=1.59e+11 * * Assumes access to binary MCPL datasets in . named [sector][beamline].mcpl.gz, i.e. W8.mcpl.gz. * @@ -30,21 +36,23 @@ * Yheight: [m] Defines the moderator height. Valid values are 0.03 m and 0.06 m * delta: [m] Parameter that allows to scan "collimator" position * thres: [ ] Weight-threshold, neutrons with higher weight are absorbed -* repeat: [ ] Number of MCPL file repetitions +* repeat: [ ] Number of MCPL file repetitions (only valid with use of MCPL_input, not with MCPL_input_once * v_smear: [1] When repeating events, make a Gaussian MC choice within v_smear*V around particle velocity V * pos_smear: [m] When repeating events, make a flat MC choice of position within pos_smear around particle starting position * dir_smear: [deg] When repeating events, make a Gaussian MC choice of direction within dir_smear around particle direction +* allmons: [1] Flag to include also debugging-monitor output +* filter: [1] Use a 'preconditioned' MCPL file that removes non-transportable neutrons (energy, position, emmision coordinates etc.) Generate such inputs using ESS_MCPL_generate_filtered. * * %L * * Benchmarking website available at http://ess_butterfly.mcstas.org * %E *******************************************************************************/ -DEFINE INSTRUMENT ESS_butterfly_MCPL_test(string sector="N", - int beamline=1,Lmin=0.2,Lmax=20,c_performance=1,t_performance=1, - int index=8,dist=0,cold=1,Yheight=0.03,delta=0,thres=0.003,filter=0,repeat=1, - v_smear=0.1,pos_smear=0.01,dir_smear=0.01) - +DEFINE INSTRUMENT ESS_butterfly_MCPL_test(string sector="S", + int beamline=2,Lmin=0.2,Lmax=20,c_performance=1,t_performance=1, + int index=0,dist=2,cold=10.5,Yheight=0.03,delta=0,thres=4e8, filter=1, repeat=1, + v_smear=0.1,pos_smear=0.01,dir_smear=0.01,int allmons=0) +SHELL " curl -O https://public.esss.dk/users/willend/MCPL/S2_filtered.mcpl.gz " DECLARE %{ double calcAlpha(double length, double radius) { @@ -73,7 +81,8 @@ DECLARE %{ double TCollmax; #pragma acc declare create(TCollmin,TCollmax) double EminTh=20, EmaxTh=100, EminC=0, EmaxC=20; - #pragma acc declare create(EminTh,EmaxTh,EminC,EmaxC) + double weightmultiplier; + #pragma acc declare create(EminTh,EmaxTh,EminC,EmaxC,weightmultiplier) /* 10 beamlines in sector N and E - plus one location added for drawing */ double iBeamlinesN[] = { 30.0, 36.0, 42.0, 48.0, 54.0, 60.0, 66.0, 72.0, 78.0, 84.0, 90.0}; double iBeamlinesE[] = {-30.0, -36.0, -42.0, -48.0, -54.0, -60.0, -66.0, -72.0, -78.0, -84.0, -90.0}; @@ -84,6 +93,8 @@ DECLARE %{ double ANGLE; double DeltaX,DeltaZ; char MCPLfile[128]; + char FILTERfile[128]; + int Divisions; %} USERVARS %{ @@ -91,8 +102,8 @@ USERVARS %{ double SrcX; double SrcY; double SrcZ; - double Emin; - double Emax; + double E_min; + double E_max; double Eneutron; double T0; double L0; @@ -135,12 +146,15 @@ INITIALIZE DeltaX=-0.0585; DeltaZ=-0.0925; } ANGLE=iBeamlines[beamline-1]-90; - if (filter==0) + if (filter==0) { sprintf(MCPLfile,"%s%i.mcpl.gz",sector,beamline); - else + printf("MCPLfile is %s\n",MCPLfile); + } else { sprintf(MCPLfile,"%s%i_filtered.mcpl.gz",sector,beamline); - printf("MCPLfile is %s\n",MCPLfile); - #pragma acc update device(EminTh,EmaxTh,EminC,EmaxC) + printf("MCPLfile is %s (FILTERED)\n",FILTERfile); + } + weightmultiplier=1.56e16/1e5; /* (ESS protons/s) / MCNP nps */ + #pragma acc update device(EminTh,EmaxTh,EminC,EmaxC,weightmultiplier) %} TRACE @@ -158,59 +172,66 @@ COMPONENT vinROT1 = Arm() AT(0,0,0) RELATIVE PREVIOUS ROTATED (-90,0,0) RELATIVE PREVIOUS - COMPONENT vin = MCPL_input(filename=MCPLfile,verbose=1,repeat_count=repeat,v_smear=v_smear,pos_smear=pos_smear,dir_smear=dir_smear) -AT(0,0,0) RELATIVE PREVIOUS -EXTEND %{ - SCATTER; - p*=1.56e16; - p/=1e5; - z=z-0.137; -%} +/* Uncomment this component to allow 'repeating' the MCPL input file */ +/* COMPONENT vin = MCPL_input(filename=MCPLfile,verbose=1,repeat_count=repeat,v_smear=v_smear,pos_smear=pos_smear,dir_smear=dir_smear) */ +/* AT(0,0,-0.137) RELATIVE PREVIOUS */ + +/* Comment this component to allow 'repeating' the MCPL input file */ +COMPONENT vin = MCPL_input_once(filename=MCPLfile,verbose=1) +AT(0,0,-0.137) RELATIVE PREVIOUS -COMPONENT Sphere1 = PSD_monitor_4PI(filename="nonrotated", radius=2.2,restore_neutron=1) -AT (0,0,0) RELATIVE PREVIOUS +COMPONENT Sphere1 = PSD_monitor_4PI(filename="nonrotated", radius=2.2,restore_neutron=1, nowritefile=!allmons) +WHEN (allmons) AT (0,0,0) RELATIVE PREVIOUS -/* Focusing for this use of the source is a little unphysical: 1x1cm @ 1m ~ 1e-4 steradian. To be useful in a "proper" instrument, you should of course illuminate your beamport fully!*/ COMPONENT Source = ESS_butterfly(sector=sector,beamline=beamline,Lmin=Lmin,Lmax=Lmax,c_performance=c_performance,t_performance=t_performance,dist=dist,target_index=index,cold_frac=cold, yheight=Yheight, focus_xw=0.12, focus_yh=0.12) WHEN (0==1) AT (DeltaX,0,DeltaZ) ABSOLUTE ROTATED (0, ANGLE, 0) ABSOLUTE -COMPONENT Sphere0 = PSD_monitor_4PI(filename="rotated", radius=2.2,restore_neutron=1) -AT (0,0,0) RELATIVE Source +COMPONENT Sphere0 = PSD_monitor_4PI(filename="rotated", radius=2.2,restore_neutron=1, nowritefile=!allmons) +WHEN (allmons) AT (0,0,0) RELATIVE Source COMPONENT Focus_cut=Shape(xwidth=0.01,yheight=0.01) - AT(0,0,2) RELATIVE Source +AT(0,0,2) RELATIVE Source EXTEND %{ + Divisions=1; double xtmp,ytmp,ztmp,vxtmp,vytmp,vztmp; xtmp=x;ytmp=y;ztmp=z; vxtmp=vx;vytmp=vy;vztmp=vz; ALLOW_BACKPROP; PROP_Z0; SCATTER; - + if (!INSTRUMENT_GETPAR(filter)) { if (fabs(x)>0.06 || fabs(y)>0.06) { ABSORB; - } else { - x=xtmp;y=ytmp;z=ztmp; - vx=vxtmp;vy=vytmp;vz=vztmp; + } + } + x=xtmp;y=ytmp;z=ztmp; + vx=vxtmp;vy=vytmp;vz=vztmp; + p*=weightmultiplier; + if (INSTRUMENT_GETPAR(thres)) { + if (p>INSTRUMENT_GETPAR(thres)) { + Divisions=10; + } } %} -COMPONENT BackTrace = Shape(xwidth=0.3,yheight=0.3) - AT (0,0,0.08) RELATIVE Source +COMPONENT ATSource=PSD_monitor(xwidth=0.31,yheight=0.31, nowritefile=!allmons) + WHEN(allmons && filter) AT (0,0,0.08) RELATIVE Source + +SPLIT Divisions COMPONENT BackTrace = Shape(xwidth=0.3,yheight=0.3) +AT (0,0,0.08) RELATIVE Source EXTEND %{ /* Propagate back to a small rectangle in front of moderators */ - + SCATTER; + double myL = (2*PI/V2K)/sqrt(vx*vx + vy*vy + vz*vz); ALLOW_BACKPROP; PROP_Z0; SCATTER; - /* Remove neutrons that are not from around the moderators */ if (fabs(x)>0.12 || fabs(y)>0.03) { ABSORB; } - double myL = (2*PI/V2K)/sqrt(vx*vx + vy*vy + vz*vz); if ( myL < INSTRUMENT_GETPAR(Lmin) || myL > INSTRUMENT_GETPAR(Lmax) ) { ABSORB; } @@ -222,16 +243,41 @@ EXTEND %{ SrcX=x;SrcY=y;SrcZ=z; Eneutron=VS2E*(vx*vx + vy*vy + vz*vz); if (Eneutron>EminTh) { - Emin=EminC;Emax=EmaxC; + E_min=EminC;E_max=EmaxC; IsCold=0; } else { - Emin=EminTh;Emax=EmaxTh; + E_min=EminTh;E_max=EmaxTh; IsCold=1; } T0=t; L0=myL; + if (Divisions>1) { + double tmpx, tmpy, tmpz; + // Position-MC: + randvec_target_circle (&tmpx, &tmpy, &tmpz, NULL, 0, 0, 1, 0); + NORM (tmpx, tmpy, tmpz); + tmpx *= INSTRUMENT_GETPAR(pos_smear) * rand01 (); + tmpy *= INSTRUMENT_GETPAR(pos_smear) * rand01 (); + tmpz *= INSTRUMENT_GETPAR(pos_smear) * rand01 (); + x += tmpx; + y += tmpy; + z += tmpz; + + double V=sqrt(vx*vx+vy*vy+vz*vz); + V *= (1 + INSTRUMENT_GETPAR(v_smear) * randpm1 ()); + double d0,d1,d2; + NORM(d0,d1,d2); + // Direction-MC: + randvec_target_circle (&d0, &d1, &d2, NULL, d0, d1, d2, sin (INSTRUMENT_GETPAR(dir_smear) * DEG2RAD)); + NORM (d0, d1, d2); + vx=d0*V; + vy=d1*V; + vz=d2*V; + } %} +COMPONENT POSTATSource=PSD_monitor(xwidth=0.31,yheight=0.31) +AT (0,0,0) RELATIVE BackTrace /* These arms are just to ensure we get a good view of the monolith */ COMPONENT Arm1 = Arm() @@ -240,51 +286,48 @@ COMPONENT Arm1 = Arm() COMPONENT Arm2 = Arm() AT (0,0,3.5) RELATIVE Source - - -COMPONENT AutoTOFL0 = Monitor_nD(xwidth=XW, yheight=YH, user1="T0", user2="L0", options="user1 limits=[0 5e-3] bins=51, user2 limits=[0.1 20] bins=41", restore_neutron=1) +COMPONENT AutoTOFL0 = Monitor_nD(xwidth=XW, yheight=YH, user1="T0", user2="L0", options="user1 limits=[0 5e-3] bins=51, user2 limits=[0.1 20] bins=41", restore_neutron=1,filename="AutoTOFL0") AT (0, 0, 0.001) RELATIVE BackTrace -COMPONENT AutoTOF0 = Monitor_nD(xwidth=XW, yheight=YH, user1="T0", options="user1 limits=[0 5e-3] bins=51", restore_neutron=1) +COMPONENT AutoTOF0 = Monitor_nD(xwidth=XW, yheight=YH, user1="T0", options="user1 limits=[0 5e-3] bins=51", restore_neutron=1,filename="AutoTOF0") AT (0, 0, 0.001) RELATIVE PREVIOUS -COMPONENT AutoL0 = Monitor_nD(xwidth=XW, yheight=YH, user1="L0", options="user1 limits=[0.1 20] bins=41", restore_neutron=1) +COMPONENT AutoL0 = Monitor_nD(xwidth=XW, yheight=YH, user1="L0", options="user1 limits=[0.1 20] bins=41", restore_neutron=1,filename="AutoL0") AT (0, 0, 0.001) RELATIVE PREVIOUS -COMPONENT PSD0= Monitor_nD(filename="flat",xwidth=0.2,yheight=0.2,user1="SrcX",user2="SrcY",options="user1 limits=[-0.1 0.1] bins=90, user2 limits=[-0.1 0.1] bins=90,", restore_neutron=1) +COMPONENT PSD0= Monitor_nD(filename="flat",xwidth=0.2,yheight=0.2,user1="SrcX",user2="SrcY",options="user1 limits=[-0.1 0.1] bins=90, user2 limits=[-0.1 0.1] bins=90,", restore_neutron=1,filename="PSD0") AT (0,0,0.001) RELATIVE PREVIOUS -COMPONENT PSD1=Monitor_nD(filename="flatC",xwidth=0.2,yheight=0.2,user1="SrcX",user2="SrcY",options="user1 limits=[-0.1 0.1] bins=90, user2 limits=[-0.1 0.1] bins=90,", restore_neutron=1) +COMPONENT PSD1=Monitor_nD(filename="flatC",xwidth=0.2,yheight=0.2,user1="SrcX",user2="SrcY",options="user1 limits=[-0.1 0.1] bins=90, user2 limits=[-0.1 0.1] bins=90,", restore_neutron=1,filename="PSD1") WHEN (Eneutron=EminTh) AT (0,0,0.001) RELATIVE PREVIOUS -COMPONENT MonND1 = Monitor_nD(xwidth=XW, yheight=YH, user1="SrcX", username1="Horizontal position / [m]", options=options1, restore_neutron=1) +COMPONENT MonND1 = Monitor_nD(xwidth=XW, yheight=YH, user1="SrcX", username1="Horizontal position / [m]", options=options1, restore_neutron=1,filename="MonND1") AT (0, 0, 1) RELATIVE Source -COMPONENT CWidth = Monitor_nD(xwidth=XW, yheight=YH, user1="SrcX", username1="Horizontal position / [m]", options=options1, restore_neutron=1) +COMPONENT CWidth = Monitor_nD(xwidth=XW, yheight=YH, user1="SrcX", username1="Horizontal position / [m]", options=options1, restore_neutron=1,filename="CWidth") WHEN(Eneutron<=EmaxC && Eneutron>=EminC) AT (0, 0, 1) RELATIVE Source -COMPONENT TWidth = Monitor_nD(xwidth=XW, yheight=YH, user1="SrcX", username1="Horizontal position / [m]", options=options1, restore_neutron=1) +COMPONENT TWidth = Monitor_nD(xwidth=XW, yheight=YH, user1="SrcX", username1="Horizontal position / [m]", options=options1, restore_neutron=1,filename="TWidth") WHEN(Eneutron<=EmaxTh && Eneutron>=EminTh) AT (0, 0, 1) RELATIVE Source -COMPONENT MonND2 = Monitor_nD(xwidth=XW, yheight=YH, user1="SrcY", username1="Vertical position COLD / [m]", options=options4, restore_neutron=1) +COMPONENT MonND2 = Monitor_nD(xwidth=XW, yheight=YH, user1="SrcY", username1="Vertical position COLD / [m]", options=options4, restore_neutron=1,filename="MonND2") WHEN(IsCold) AT (0, 0, 1) RELATIVE Source -COMPONENT MonND2_2 = Monitor_nD(xwidth=XW, yheight=YH, user1="SrcY", username1="Vertical position THERMAL/ [m]", options=options4, restore_neutron=1) +COMPONENT MonND2_2 = Monitor_nD(xwidth=XW, yheight=YH, user1="SrcY", username1="Vertical position THERMAL/ [m]", options=options4, restore_neutron=1,filename="MonND2_2") WHEN(!IsCold) AT (0, 0, 1) RELATIVE Source -COMPONENT MonND3 = Monitor_nD(xwidth=XW, yheight=YH, user1="SrcX", username1="Horizontal position / [m]", user2="SrcY",username2="Vertical position / [m]", options=options2, restore_neutron=1) +COMPONENT MonND3 = Monitor_nD(xwidth=XW, yheight=YH, user1="SrcX", username1="Horizontal position / [m]", user2="SrcY",username2="Vertical position / [m]", options=options2, restore_neutron=1,filename="MonND3") AT (0, 0, 1) RELATIVE Source -COMPONENT MonND4 = Monitor_nD(xwidth=XW, yheight=YH, user1="SrcX", username1="Emission position / [m]", user2="SrcZ", username2="Z-component of position / [m]", options="user1 bins=201 limits=[-0.3,0.3], user2 bins=201 limits=[-0.3,0.3]", restore_neutron=1) +COMPONENT MonND4 = Monitor_nD(xwidth=XW, yheight=YH, user1="SrcX", username1="Emission position / [m]", user2="SrcZ", username2="Z-component of position / [m]", options="user1 bins=201 limits=[-0.3,0.3], user2 bins=201 limits=[-0.3,0.3]", restore_neutron=1,filename="MonND4") AT (0, 0, 1) RELATIVE Source - -COMPONENT AutoTOFL = Monitor_nD(xwidth=0.1, yheight=0.1, options="tof limits=[0 15e-3] bins=51, lambda limits=[0.1 20] bins=41", restore_neutron=1) +COMPONENT AutoTOFL = Monitor_nD(xwidth=0.1, yheight=0.1, options="tof limits=[0 15e-3] bins=51, lambda limits=[0.1 20] bins=41", restore_neutron=1,filename="AutoTOFL") AT (0, 0, 2) RELATIVE Source @@ -293,8 +336,8 @@ COMPONENT AutoTOFL = Monitor_nD(xwidth=0.1, yheight=0.1, options="tof limits=[0 COMPONENT BrillmonCOLD = Brilliance_monitor( nlam = 101, nt = 101, filename = "brillCOLD", t_0 = -1000, t_1 =4e4, lambda_0 = lambdamin, lambda_1 = lambdamax, - Freq =14, toflambda=1 ,tofcuts=0, srcarea=(100*0.072*100*Yheight), restore_neutron=1,source_dist=2,xwidth=0.1,yheight=0.1) -WHEN(IsCold) AT (0, 0, 2) RELATIVE Source + Freq =14, toflambda=1 ,tofcuts=0, srcarea=(100*0.072*100*Yheight), restore_neutron=1,source_dist=2,xwidth=0.1,yheight=0.1,nowritefile=!allmons) +WHEN(allmons && IsCold) AT (0, 0, 2) RELATIVE Source /* Measures "collimated" brilliance of the cold source over fixed 6 cm wide area x central part vertically. */ /* Used for calibration of performance wrt. MCNP BF1 output, see http://ess_butterfly.mcstas.org */ @@ -308,8 +351,8 @@ COMPONENT BrillmonCOLD_COLL = Brilliance_monitor( COMPONENT BrillmonTHRM = Brilliance_monitor( nlam = 101, nt = 101, filename = "brillTHRM", t_0 = -1000, t_1 =4e4, lambda_0 = lambdamin, lambda_1 = lambdamax, - Freq =14, toflambda=1,tofcuts=0, srcarea=(100*0.108*100*Yheight), restore_neutron=1,source_dist=2,xwidth=0.1,yheight=0.1) - WHEN (!IsCold) AT (0, 0, 2) RELATIVE Source + Freq =14, toflambda=1,tofcuts=0, srcarea=(100*0.108*100*Yheight), restore_neutron=1,source_dist=2,xwidth=0.1,yheight=0.1,nowritefile=!allmons) + WHEN (allmons && !IsCold) AT (0, 0, 2) RELATIVE Source /* Measures "collimated" brilliance of the thermal source over fixed 6 cm wide area (or smaller at beamlines no. 1,2) x central part vertically. */ /* Used for calibration of performance wrt. MCNP BF1 output, see http://ess_butterfly.mcstas.org */ @@ -325,12 +368,20 @@ COMPONENT PSD0x= Monitor_nD(filename="flat_x",xwidth=0.1,yheight=0.1,user1="SrcX COMPONENT PSD1x=Monitor_nD(filename="flatC_x",xwidth=0.1,yheight=0.1,user1="SrcX",user2="SrcY",options="user1 limits=[-0.08 0.08] bins=90, user2 limits=[-0.08 0.08] bins=90,", restore_neutron=1) WHEN (Eneutron=EminTh) AT (0,0,0.001) RELATIVE PREVIOUS +COMPONENT Flux_incoming=Monitor_nD(xwidth=0.05,yheight=0.1,filename="Flux_incoming",options="x, y pr cm2", bins=100) + AT (0, 0, 2) RELATIVE Source + +COMPONENT GuideStraight=Guide(w1 = 0.05, h1 = 0.1, l = 10) + AT (0,0,0) RELATIVE PREVIOUS + +COMPONENT Flux_mid=Monitor_nD(xwidth=0.05,yheight=0.1,filename="Flux_mid",options="x, y pr cm2", bins=100) + AT (0, 0, 10) RELATIVE GuideStraight COMPONENT GuideR = Guide_curved( - w1 = 0.1, h1 = 0.1, l = 50, curvature=3000) + w1 = 0.05, h1 = 0.1, l = 50, curvature=3000) AT (0,0,0.001) RELATIVE PREVIOUS COMPONENT RArm=Arm() @@ -342,9 +393,32 @@ COMPONENT GuideR = Guide_curved( ROTATED (0, calcAlpha(50,3000),0) RELATIVE GuideR COMPONENT Monitor2_xy1 = Monitor_nD( - options = "x limits=[-0.06 0.06] bins=51, y limits=[-0.06 0.06] bins=51,", xwidth = 0.12, yheight = 0.12) + options = "x limits=[-0.06 0.06] bins=51, y limits=[-0.06 0.06] bins=51,", xwidth = 0.12, yheight = 0.12,filename="Monitor2_xy1") AT (0, 0, 0.05) RELATIVE RArm +COMPONENT Monitor_X = Monitor_nD( + options = "x limits=[-0.06 0.06] bins=51", xwidth = 0.12, yheight = 0.12,filename="Monitor_X") + AT (0, 0, 0) RELATIVE Monitor2_xy1 + +COMPONENT Monitor_Y = Monitor_nD( + options = "y limits=[-0.06 0.06] bins=51", xwidth = 0.12, yheight = 0.12,filename="Monitor_Y") + AT (0, 0, 0) RELATIVE Monitor2_xy1 + +COMPONENT Monitor_divH = Monitor_nD( + options = "hdiv limits=[-1.5 1.5] bins=51", xwidth = 0.12, yheight = 0.12,filename="Monitor_divH") + AT (0, 0, 0) RELATIVE Monitor2_xy1 + +COMPONENT Monitor_divV = Monitor_nD( + options = "vdiv limits=[-1.5 1.5] bins=51", xwidth = 0.12, yheight = 0.12,filename="Monitor_divV") + AT (0, 0, 0) RELATIVE Monitor2_xy1 + +COMPONENT Monitor_t = Monitor_nD( + options = "t limits=[0 350e-3] bins=351", xwidth = 0.12, yheight = 0.12,filename="Monitor_t") + AT (0, 0, 0) RELATIVE Monitor2_xy1 + +COMPONENT Flux_end=Monitor_nD(xwidth=0.05,yheight=0.1,filename="Flux_end",options="x, y pr cm2", bins=100) + AT (0, 0, 0) RELATIVE Monitor2_xy1 + /* /\* Uncomment these helper-arms to view "full" monolith *\/ */ COMPONENT DummyArm1 = Arm() diff --git a/mcstas-comps/examples/ESS/ESS_butterfly_MCPL_test/README.md b/mcstas-comps/examples/ESS/ESS_butterfly_MCPL_test/README.md index 93c96c4a8..288bc128c 100644 --- a/mcstas-comps/examples/ESS/ESS_butterfly_MCPL_test/README.md +++ b/mcstas-comps/examples/ESS/ESS_butterfly_MCPL_test/README.md @@ -1,6 +1,6 @@ # The `ESS_butterfly_MCPL_test` Instrument -*McStas: Test instrument for the updated BF1 butterfly moderator design using MCPL input* +*McStas: Test instrument for the ESS source using direct MCPL input. Includes a simple guide.* ## Identification @@ -12,12 +12,14 @@ ## Description ```text -Test instrument for the updated BF1 butterfly moderator design using MCPL input files. - -The below example gives a 50-50 (statistics-wise) cold/thermal beam at beamline N10. -Example: sector=N beamline=10 cold=0.5 - -Assumes access to binary MCPL datasets in . named [sector][beamline].mcpl.gz, i.e. W8.mcpl.gz. +Test instrument for the updated BF1 butterfly moderator design using MCPL input files. The instrument includes a simple guide with 10 m straight flight-path followed by a 50 m curved guide. + +The below example gives a 50-50 (statistics-wise) cold/thermal beam at beamline S2 using the base S2 MCPL file from +https://public.esss.dk/users/willend/MCPL/S2.mcpl.gz +Example: sector=S beamline=2 Detector: Monitor2_xy1_I=1.59e+11 +and a 'preconditioned' version of this the file from https://public.esss.dk/users/willend/MCPL/S2_filtered.mcpl.gz +(This file is fetched via a curl command at code-generation time) +Example: sector=S beamline=2 filter=1 Detector: Monitor2_xy1_I=1.59e+11 ``` ## Input parameters @@ -26,23 +28,24 @@ Parameters in **boldface** are required; the others are optional. | Name | Unit | Description | Default | |------|------|-------------|---------| -| sector | str | Defines the 'sector' of your instrument position. Valid values are "N","S","E" and "W" | "N" | -| beamline | 1 | Defines the 'beamline number' of your instrument position. Valid values are 1..10 or 1..11 depending on sector | 1 | +| sector | str | Defines the 'sector' of your instrument position. Valid values are "N","S","E" and "W" | "S" | +| beamline | 1 | Defines the 'beamline number' of your instrument position. Valid values are 1..10 or 1..11 depending on sector | 2 | | Lmin | AA | Minimum wavelength simulated | 0.2 | | Lmax | AA | Maximum wavelength simulated | 20 | | c_performance | 1 | Cold brilliance scalar performance multiplicator c_performance > 0 | 1 | | t_performance | 1 | Thermal brilliance scalar performance multiplicator t_performance > 0 | 1 | -| index | 1 | Target index for source focusing. Defaults to illuminate the "cold collimated" brilliance monitor, thereby suppressing "dist" | 8 | -| dist | m | Distance from origin to focusing rectangle; at (0,0,dist) - alternatively use target_index | 0 | -| cold | 1 | Defines the statistical fraction of events emitted from the cold part of the moderator | 1 | +| index | 1 | Target index for source focusing. Defaults to illuminate the "cold collimated" brilliance monitor, thereby suppressing "dist" | 0 | +| dist | m | Distance from origin to focusing rectangle; at (0,0,dist) - alternatively use target_index | 2 | +| cold | 1 | Defines the statistical fraction of events emitted from the cold part of the moderator | 10.5 | | Yheight | m | Defines the moderator height. Valid values are 0.03 m and 0.06 m | 0.03 | | delta | m | Parameter that allows to scan "collimator" position | 0 | -| thres | | Weight-threshold, neutrons with higher weight are absorbed | 0.003 | -| filter | | | 0 | -| repeat | | Number of MCPL file repetitions | 1 | +| thres | | Weight-threshold, neutrons with higher weight are absorbed | 4e8 | +| filter | 1 | Use a 'preconditioned' MCPL file that removes non-transportable neutrons (energy, position, emmision coordinates etc.) Generate such inputs using ESS_MCPL_generate_filtered. | 1 | +| repeat | | Number of MCPL file repetitions (only valid with use of MCPL_input, not with MCPL_input_once | 1 | | v_smear | 1 | When repeating events, make a Gaussian MC choice within v_smear*V around particle velocity V | 0.1 | | pos_smear | m | When repeating events, make a flat MC choice of position within pos_smear around particle starting position | 0.01 | | dir_smear | deg | When repeating events, make a Gaussian MC choice of direction within dir_smear around particle direction | 0.01 | +| allmons | 1 | Flag to include also debugging-monitor output | 0 | ## Links