/******************************************************************************* * Instrument: ESS_butterfly_Guide_curved_test * * %I * Written by: Peter Willendrup * Date: 2016-09-26 * Origin: ESS * %INSTRUMENT_SITE: ESS * * Test instrument for the updated BF1 butterfly moderator design including a simple guide. * * %D * 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 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 * 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 * 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="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 %{ 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; #pragma acc declare create(EminTh,EmaxTh,EminC,EmaxC) /* 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; %} USERVARS %{ int IsCold; double SrcX; double SrcY; double SrcZ; double E_min; double E_max; double Eneutron; double T0; double L0; %} 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; #pragma acc update device(EminTh,EmaxTh,EminC,EmaxC) %} 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 COMPONENT vinROT1 = Arm() AT(0,0,0) RELATIVE PREVIOUS ROTATED (-90,0,0) RELATIVE PREVIOUS COMPONENT vin = Arm() 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, 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, 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 EXTEND %{ double myL = (2*PI/V2K)/sqrt(vx*vx + vy*vy + vz*vz); /* Measure location and energy for later use */ SrcX=x;SrcY=y;SrcZ=z; Eneutron=VS2E*(vx*vx + vy*vy + vz*vz); if (Eneutron>EminTh) { E_min=EminC;E_max=EmaxC; IsCold=0; } else { 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() AT (0,0,2) RELATIVE Source 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,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,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,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,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,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,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,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,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,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,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,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,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,filename="AutoTOFL") AT (0, 0, 2) RELATIVE Source /* Measures brilliance of the "full" cold source */ 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,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 */ COMPONENT BrillmonCOLD_COLL = Brilliance_monitor( nlam = 101, nt = 101, filename = "brillCOLD_COLL", t_0 = 0, t_1 = 4e4, lambda_0 = lambdamin, lambda_1 = lambdamax, Freq =14, toflambda=1,tofcuts=0, srcarea=(100*0.06*100*2*Yheight/2.5), restore_neutron=1,source_dist=2,xwidth=0.1,yheight=0.1) WHEN(fabs(SrcY) (0.011+delta)) AT (0, 0, 2) RELATIVE Source /* Measures brilliance of the "full" thermal source */ 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,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 */ COMPONENT BrillmonTHRM_COLL = Brilliance_monitor( nlam = 101, nt = 101, filename = "brillTHRM_COLL", t_0 = -1000, t_1 =4e4, lambda_0 = lambdamin, lambda_1 = lambdamax, Freq =14, toflambda=1,tofcuts=0, srcarea=(100*0.06*100*2*Yheight/2.5), restore_neutron=1,source_dist=2,xwidth=0.1,yheight=0.1) WHEN (fabs(SrcY)(TCollmin+delta) && -SrcX<(TCollmax+delta)) AT (0, 0, 2) RELATIVE Source COMPONENT PSD0x= Monitor_nD(filename="flat_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) AT (0,0,0.001) RELATIVE PREVIOUS 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.05, h1 = 0.1, l = 50, curvature=3000) AT (0,0,0.001) RELATIVE PREVIOUS COMPONENT RArm=Arm() AT (0.416657,0,50.002) RELATIVE GuideR ROTATED (0, calcAlpha(50,3000),0) RELATIVE GuideR COMPONENT RArm2=Arm() AT (calcX(50, 3000),0,calcZ(50, 3000)) RELATIVE GuideR 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,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() AT (6,0,6) ABSOLUTE COMPONENT DummyArm2 = Arm() AT (-6,0,6) ABSOLUTE COMPONENT DummyArm3 = Arm() AT (-6,0,-6) ABSOLUTE COMPONENT DummyArm4 = Arm() AT (6,0,-6) ABSOLUTE COMPONENT DummyArm5 = Arm() AT (6,0,6) ABSOLUTE FINALLY %{ %} END