18 Commits

Author SHA1 Message Date
glavic_a c85792ef22 Update model from DMSC version for McStas 3.x 2025-10-24 14:48:27 +02:00
glavic_a 6e0d611e73 add pycharm to gitignore 2025-10-24 14:37:09 +02:00
glavic_a b570182e5d Fix feeder not beeing fully illuminated and add VS and MF divergence monitors 2021-10-13 11:39:24 +02:00
glavic_a 41533faaa1 Add polarization parameters of SNAG mirrors and run file 2021-10-13 10:00:11 +02:00
glavic_a 1cbe0c628f Implement geometry of double V-style polarizer, straight not correct reflection, yet 2021-10-11 22:55:53 +02:00
glavic_a c2a947a923 Minor changes to run on cluster 2019-11-18 16:45:24 +01:00
glavic_a c011bdad07 First version for polarizer distance sweep 2019-08-02 19:23:48 +02:00
glavic_a f6a77d0dfe Fix issue with multiple polarizer instances reintroduced in update 2019-07-28 19:23:31 +02:00
glavic_a 7323c7fd8a Remove monitors, add geometry parameters and slurm job. Some strange behavior 2019-07-28 17:40:21 +02:00
glavic_a fb84f6fb5b Update polarizer and optimize parameters for m=5 mirror from SNAG 2019-07-28 14:50:31 +02:00
glavic_a fe44ee1edd Add updated Polarizer component 2019-07-25 15:40:51 +02:00
glavic_a 6016ffba54 Replace Selene guide top/bottom order 2019-07-25 15:40:28 +02:00
glavic_a 0cc9f0ae5d Add monitor model and update CPC1 entrance and middle geometry 2018-09-15 14:13:56 +02:00
glavic_a 00bd8b18c0 Fix reader plotting for 2D datasets 2018-09-15 14:11:19 +02:00
ll_glavic 7a49a0f22c Added selene guide position error based on foot missalignment and MCC Slurm job so sweep the parameters. 2018-08-18 11:15:15 +02:00
glavic_a a58b2d3030 Separate Selene 1/2 and middle focus 2018-07-03 08:22:36 +02:00
glavic_a 6dee3c98c1 Start monitor cleanup by removing all put the essential nD detectors. 2018-07-03 08:21:55 +02:00
glavic_a 8c10d63987 Selene guide m-value per segment 2018-06-26 17:53:46 +02:00
29 changed files with 2437 additions and 1168 deletions
+1
View File
@@ -8,6 +8,7 @@ old
*.backup
.project
.settings
.idea
*.pyc
*.pyo
hosts
+9
View File
@@ -0,0 +1,9 @@
#!/bin/bash
# compress McStas hdf5 files to increase performance and save disk space
for fi in $*
do
echo $fi
h5repack -i $fi -o $fi.c -f /entry1/data/detector_list_p_x_y_t_L_sx_sy_sz/events:GZIP=5 -l /entry1/data/detector_list_p_x_y_t_L_sx_sy_sz/events:CHUNK=3072x8
mv $fi.c $fi
done
+20 -20
View File
@@ -11,7 +11,7 @@ from numpy import *
try:
import h5py
except ImportError:
print "h5py not found, modern NeXuS format will not be readable."
print("h5py not found, modern NeXuS format will not be readable.")
try:
from IPython import display #@UnusedImport
@@ -49,7 +49,7 @@ class McSim(object):
elif os.path.exists(os.path.join(path, 'mccode.sim')):
self._init_old(os.path.join(path, 'mccode.sim'))
else:
raise IOError, "Can't locate mccode.h5 of mccode.sim file in %s"%path
raise IOError("Can't locate mccode.h5 of mccode.sim file in %s"%path)
def _init_hdf(self, path):
self.hdf=h5py.File(path, 'r')
@@ -61,11 +61,11 @@ class McSim(object):
self.data_loader=DataLoaderOld(self.info, os.path.dirname(path))
def keys(self):
return self.info['data'].keys()
return list(self.info['data'].keys())
def monitors(self):
output=[]
for val in self.info['data'].values():
for val in list(self.info['data'].values()):
xy=(val['xvar'], val['yvar'])
if not xy in output:
output.append(xy)
@@ -74,7 +74,7 @@ class McSim(object):
def plot(self, monitors=None):
graphs={}
for key, val in self.info['data'].items():
for key, val in list(self.info['data'].items()):
xy=(val['xvar'], val['yvar'])
if monitors is not None and xy!=monitors:
continue
@@ -102,12 +102,12 @@ class McSim(object):
def __getitem__(self, item):
if item in self._data:
return self._data[item]
elif item in self.keys():
elif item in list(self.keys()):
data=self.data_loader.load_item(item)
self._data[item]=data
return data
else:
raise KeyError, "Can't find dataset %s"%item
raise KeyError("Can't find dataset %s"%item)
class HeaderFile(object):
'''
@@ -124,7 +124,7 @@ class HeaderFile(object):
data=open(path, 'r').read()
data_lines=data.splitlines()
if not data.startswith('McStas simulation description file'):
raise IOError, 'Not a valid McStas description file.'
raise IOError('Not a valid McStas description file.')
self._data['start_time']=data_lines[1].split(':', 1)[1].strip()
self._data['program_name']=data_lines[2].split(':', 1)[1].strip()
self._data['data']=self.get_data(data)
@@ -167,13 +167,13 @@ class DataLoaderOld(object):
fname=os.path.join(self.root, item_info['filename'])
x_col=item_info['xvar']
y_col=item_info['yvar']
if x_col=='Li' and y_col=='p': # Detector_nD
if x_col in ['Li', 'List'] and y_col=='p': # Detector_nD
cols=item_info['variables'].split()
data=loadtxt(fname, dtype={'names': cols, 'formats': ['f4']*len(cols)})
return TofData(data, item_info)
else:
raw=loadtxt(fname)
data=raw[:len(raw)/3]
data=raw[:len(raw)//3]
return Dataset(data, item_info)
def load_item_1d(self, item):
@@ -192,11 +192,11 @@ class DataLoaderHDF(object):
self.hdf=hdf['entry1']
self.info={}
self.info['data']={}
for item in self.hdf['data'].keys():
for item in list(self.hdf['data'].keys()):
node=self.hdf['data/'+item]
info={}
for key, value in node.attrs.items():
info[key.strip()]=value.strip()
for key, value in list(node.attrs.items()):
info[key.strip()]=value.strip().decode('ascii')
if info['filename'][-4:]=='.dat' or '_list.' in info['filename']:
self.info['data'][info['component']]=info
else:
@@ -211,7 +211,7 @@ class DataLoaderHDF(object):
node=self.hdf[item_info['datapath']]
x_col=item_info['xvar']
y_col=item_info['yvar']
if x_col=='Li' and y_col=='p': # Detector_nD
if x_col in ['Li', 'List'] and y_col=='p': # Detector_nD
cols=item_info['variables'].split()
evds=node['events']
if len(evds)<=MAX_EVTS_BATCH:
@@ -253,7 +253,7 @@ class Dataset1D(object):
import pylab
ax=pylab.gca()
limits=map(float, self.info['xlimits'].split())
limits=list(map(float, self.info['xlimits'].split()))
x=linspace(limits[0], limits[1], len(self.data))
ax.errorbar(x, self.data, yerr=self.errors)
@@ -291,12 +291,12 @@ class Dataset(object):
import pylab
ax=pylab.gca()
limits=map(float, self.info['xylimits'].split())
limits=list(map(float, self.info['xylimits'].split()))
if log:
ax.imshow(self.data, extent=limits, aspect='auto', norm=LogNorm())
ax.imshow(self.data, origin='lower', extent=limits, aspect='auto', norm=LogNorm())
else:
ax.imshow(self.data, extent=limits, aspect='auto')
ax.imshow(self.data, origin='lower', extent=limits, aspect='auto')
ax.set_xlabel(self.info['xlabel'])
ax.set_ylabel(self.info['ylabel'])
ax.set_title(self.info['component'])
@@ -342,7 +342,7 @@ class TofData(Dataset):
if fltr is None:
I, x=histogram(columns[col], bins=bins, weights=w)
else:
if isinstance(fltr, basestring):
if isinstance(fltr, str):
fltr=eval(fltr, globals(), columns)
I, x=histogram(columns[col][fltr], bins=bins, weights=w[fltr])
return x, I
@@ -367,7 +367,7 @@ class TofData(Dataset):
I, y, x=histogram2d(columns[ycol], columns[xcol],
bins=bins, weights=self.data['p'])
else:
if isinstance(fltr, basestring):
if isinstance(fltr, str):
fltr=eval(fltr, globals(), columns)
I, y, x=histogram2d(columns[ycol][fltr], columns[xcol][fltr],
bins=bins, weights=columns['p'][fltr])
+97 -175
View File
@@ -53,8 +53,10 @@ DEFINE INSTRUMENT ESS_reflectometer_Estia
int operationmode = 0, double over_illumination = 0.0, double theta_resolution = 0.04,
double lambda_min = 3.75, double lambda_start = 3.0, double lambda_end = 12.0,
int enable_chopper = 0, int enable_gravity=0, int enable_windows=1,
int enable_polarizer = 0, int enable_analyzer = 0,
double source_power = 2, double frame_usage = 0.98
int enable_polarizer = 0, int enable_analyzer = 0,
double pol1_start=0.3, double pol1_angle=1.66, double pol2_start=0.3, double pol2_angle=1.66,
double source_power = 2,
double selene1_foot1y =0.0, double selene1_foot2y = 0.0
)
DECLARE
@@ -70,8 +72,12 @@ double iscs_z=0.0547095; // ISCS in McStas coordinates using mod-view-opt.instr
double iscs_rot_y=0.4; // ISCS is at TCS-35.6 degree, McStas at -36 degree
double iscs_rot_x=0.7; // downward tilt of Estia axis
// Selene 1 (Will follow later with detailed parameters)
// Selene 1 geometry parameters (optics parameters in Estia_selene1.instr)
double selene1_foot1 = 4.20; // distance of first foot to VS focus
double selene1_foot2 = 7.00; // distance of second foot to VS focus
double selene1_center;
double selene1_shift;
double selene1_rot;
// Selene 2
@@ -97,30 +103,35 @@ double slit_distance = 1.775; //m, distance slit to sample
double selene_theta = 1.25;
double velocity_max ; // m/s neutron velocity of lambda_min
/* flags */
int PP_small = 1 ; // plot all PSD with small area
int PP_large = 1 ; // plot all PSD with large area
double analyzer_max_height = 0.01; // Maximum sample height to be covered by te polarization analyzers
double analyzer_flipper_start = 0.4; // Distance from sample to analyzer flipper
double analyzer1_start = 0.65; // Distance to sample to start the first analyzer
double analyzer1_length = 1.3; // length of first analyzer
double analyzer2_start = 0.7; // Distance to sample to start the first analyzer
double analyzer2_length = 1.6; // length of first analyzer
double Theta1_analyzer1, Theta1_analyzer2, Theta2_analyzer1, Theta2_analyzer2, dist_ana_vfocus; // quantities calculated out of values above and 1.5 degree covered divergence
%}
USERVARS %{
int p_int; // a flag that gets incremented if a polarizer mirror scatters
%}
INITIALIZE
%{
// chopper coupling together
velocity_max= 3.956034E3 / lambda_min; // h/m_n over lambda - ((3.9..e-7m^2/s))/(1e-10m)
chopper_phase=360.0*(chopper_pos*chopper_freq/velocity_max+(pulse_zero-opening_time)*chopper_freq);
chopper_open =360.0*(chopper_pos/(total_length+detector_arm))*frame_usage;
chopper_open = 98.0;
/* print out some calculated parameter for checking purposes */
printf(" Chopper phase = %.1f deg\n", chopper_phase);
printf(" Chopper open = %.1f deg\n", chopper_open);
//printf(" Chopper open = %.1f deg\n", chopper_open);
dist_ana_vfocus=analyzer_max_height/2.0/tan(1.5*PI/360.0); // the virtual focus point infront of the actual sample focus where the beams furthest out meet
@@ -130,7 +141,11 @@ Theta2_analyzer1=atan((dist_ana_vfocus+analyzer1_start+analyzer1_length)/dist_an
Theta1_analyzer2=atan((dist_ana_vfocus+analyzer2_start)/dist_ana_vfocus*analyzer_max_height/2.0/analyzer2_start)*180.0/PI;
Theta2_analyzer2=atan((dist_ana_vfocus+analyzer2_start+analyzer2_length)/dist_ana_vfocus*analyzer_max_height/2.0/(analyzer2_start+analyzer2_length))*180.0/PI;
selene1_center=(selene1_foot1+selene1_foot2)/2.;
selene1_shift=(selene1_foot1y+selene1_foot2y)/1000.0;
selene1_rot=atan((selene1_foot2y-selene1_foot1y)/(selene1_foot2-selene1_foot1)/1000.0)*180.0/PI;
//printf(" Selene 1 rotation = %.4f deg\n", selene1_rot);
//printf(" Selene 1 shift = %.1f mm\n", selene1_shift*1e3);
%}
TRACE
@@ -157,10 +172,12 @@ COMPONENT ISCS = Arm() // rotate around y-axis (slight downward tilt)
COMPONENT arm_feeder = Arm()
AT (0, 0, 0) RELATIVE ISCS
ROTATED (0, 0, 0) RELATIVE ISCS
COMPONENT arm_selene1_center = Arm()
AT (selene1_shift, 0, 2*NBOA_c+selene1_center) RELATIVE ISCS
ROTATED (0, selene1_rot, 0) RELATIVE ISCS
// Axes starting at focus of feeder parallel to c-axis of Selene guides
COMPONENT arm_selene1 = Arm()
AT (0, 0, 2*NBOA_c) RELATIVE ISCS
ROTATED (0, 0, 0) RELATIVE ISCS
AT (0, 0, -selene1_center) RELATIVE arm_selene1_center
COMPONENT arm_selene2 = Arm()
AT (0, 0, 2*NBOA_c+2*selene_c) RELATIVE ISCS
ROTATED (0, 0, 0) RELATIVE ISCS
@@ -194,20 +211,11 @@ COMPONENT moderator = ESS_butterfly(
sector = "E", beamline = 2, yheight = 0.03, cold_frac = 0.9,
focus_xw = E02_01_01_Cu_in_xmax-E02_01_01_Cu_in_xmin+0.002,
focus_yh = E02_01_01_Cu_in_yheight+0.002,
target_index=4, Lmin = lambda_start, Lmax = lambda_end,
target_index=3, Lmin = lambda_start, Lmax = lambda_end,
n_pulses = 1+enable_chopper, acc_power=source_power)
AT (0, 0, 0) RELATIVE origin
// monitor used for brilliance transfer calculations
COMPONENT tof_source = Monitor_nD(
filename = "tof_source",
options = "xdiv limits=[-0.75 0.75] bins=15 ydiv limits=[-2.0 2.0] bins=40 time limits=[0 0.6] bins=600 lambda limits=[0 35] bins=350",
xwidth=0.01, yheight = 0.01)
WHEN sample==4
AT (0, 0, 0.03) RELATIVE ISCS
ROTATED (0,-selene_theta,0) RELATIVE ISCS
/***************************************
* Geometry of neutron feeder separate *
***************************************/
@@ -216,11 +224,14 @@ COMPONENT tof_source = Monitor_nD(
/****************************************************
* Beam manipulation area around the virtual source *
****************************************************/
/* Absorber to reduce beam to needed size an for shielding purposes (FeNi in CAD model) */
COMPONENT FeNi_in = Slit(xwidth=0.028, yheight=0.078)
AT (0, 0, -0.800) RELATIVE arm_virtual_source_beam
/* Absorber to reduce beam to needed size an for shielding purposes (CPC1 in CAD model) */
COMPONENT CPC1_in = Slit(xwidth=0.0303, yheight=0.0792)
AT (0, 0, -0.890) RELATIVE arm_virtual_source_beam
COMPONENT FeNi_out = Slit(xwidth=0.013, yheight=0.038)
COMPONENT CPC1_monitor = Slit(xwidth=0.016344, yheight=0.05547)
AT (0, 0, -0.3573) RELATIVE arm_virtual_source_beam
COMPONENT CPC1_out = Slit(xwidth=0.013, yheight=0.038)
AT (0, 0, -0.220) RELATIVE arm_virtual_source_beam
@@ -231,31 +242,16 @@ COMPONENT chopper = DiskChopper(radius=chopper_diameter/2.0, yheight=0.02,
WHEN enable_chopper==1
AT (0, 0, chopper_pos-2*NBOA_c) RELATIVE arm_virtual_source_beam
// Beam monitoring before VS
COMPONENT DP_before_virtual_source = DivPos_monitor(
filename = "DP_before_virtual_source" ,
nh = 200, ndiv = 200,
xwidth = 0.2, yheight = 0.2,
maxdiv_h = 6,
restore_neutron = 1)
AT (0, 0, -0.022) RELATIVE arm_virtual_source_beam
COMPONENT PP_small_before_virtual_source = PSD_monitor(
filename = "PP_small_before_virtual_source",
nx = 200, ny = 200,
xwidth = 0.1, yheight = 0.1,
restore_neutron = 1)
WHEN PP_small ==1
AT (0, 0, -0.021) RELATIVE arm_virtual_source_beam
COMPONENT PP_large_before_virtual_source = PSD_monitor(
filename = "PP_large_before_virtual_source",
nx = 200, ny = 200,
xmin = -1.8, xmax = 1.8, ymin = -1.8, ymax = 1.8,
restore_neutron = 1)
WHEN PP_large ==1
AT (0, 0, 2*NBOA_c-0.02) RELATIVE ISCS
COMPONENT vs_divergence_h = DivPos_monitor(nb=21, ndiv=41, filename="vs_hordiv",
xmin=-0.02, xmax=0.02, ymin=-0.02, ymax=0.02, maxdiv=2.0,
restore_neutron=1)
AT (0, 0, -0.5*sample_length-0.001) RELATIVE arm_selene1
COMPONENT vs_divergence_v = DivPos_monitor(nb=21, ndiv=41, filename="vs_verdiv",
xmin=-0.02, xmax=0.02, ymin=-0.02, ymax=0.02, maxdiv=2.0,
restore_neutron=1, vertical=1)
AT (0, 0, -0.5*sample_length-0.001) RELATIVE arm_selene1
/* The actual virtual source mask, two L-shaped absorbers (first top-right) */
COMPONENT virtual_source_TR = Slit(
@@ -263,66 +259,26 @@ COMPONENT virtual_source_TR = Slit(
WHEN sample!=4
AT (-over_illumination, 0, -0.5*sample_length) RELATIVE arm_virtual_source
// window to absorb neutron not passing through heavy collimation
COMPONENT virtual_source_HC = Slit(
xmin = -0.015, xmax = 0.015,
ymin = -0.015, ymax = 0.015)
// window to cut down to defined size for test setting
COMPONENT virtual_source_HC = Slit(xwidth=sample_length, yheight=sample_height)
WHEN sample==4
AT (0, 0, 0) RELATIVE arm_virtual_source
//
/* The actual virtual source mask, two L-shaped absorbers (second bottom-left) */
COMPONENT virtual_source_BL = Slit(
xmin = -1.0, xmax = 0.0, ymin = -sample_height/2-over_illumination*5, ymax = 1.0)
WHEN sample!=4
AT (over_illumination, 0, 0.5*sample_length) RELATIVE arm_virtual_source
// Beam monitoring after VS
COMPONENT tof_virtual_source = Monitor_nD(
filename = "tof_virtual_source",
options = "xdiv limits=[-0.75 0.75] bins=15 ydiv limits=[-2.0 2.0] bins=40 time limits=[0 0.6] bins=600 lambda limits=[0 35] bins=350",
xwidth=sample_length, yheight = sample_height)
WHEN sample==4
AT (0, 0, 0.019) RELATIVE arm_virtual_source_beam
COMPONENT DL_behind_virtual_source = DivLambda_monitor(
filename = "DL_behind_virtual_source",
nL = 200, nh = 200,
xwidth = 0.02, yheight = 0.06,
maxdiv_h = 3,
Lmin = lambda_start, Lmax = lambda_end,
restore_neutron = 1)
AT (0, 0, 0.02) RELATIVE arm_virtual_source_beam
COMPONENT DP_behind_virtual_source = DivPos_monitor(
filename = "DP_behind_virtual_source" ,
nh = 200, ndiv = 200,
xwidth = 0.2, yheight = 0.2,
maxdiv_h = 3,
restore_neutron = 1)
AT (0, 0, 0.021) RELATIVE arm_virtual_source_beam
COMPONENT PP_small_behind_virtual_source = PSD_monitor(
filename = "PP_small_behind_virtual_source",
nx = 200, ny = 200,
xwidth = 0.02, yheight = 0.05,
restore_neutron = 1)
WHEN PP_small ==1
AT (0, 0, 0.022) RELATIVE arm_virtual_source_beam
COMPONENT PP_large_behind_virtual_source = PSD_monitor(
filename = "PP_large_behind_virtual_source",
nx = 200, ny = 200,
xmin = -1.8, xmax = 1.8, ymin = -1.8, ymax = 1.8,
restore_neutron = 1)
WHEN PP_large ==1
AT (0, 0, 2*NBOA_c+0.023) RELATIVE ISCS
/*************************************
* Geometry of Selene guide separate *
*************************************/
%include "Estia_selene.instr"
%include "Estia_selene1.instr"
%include "Estia_mf.instr"
%include "Estia_selene2.instr"
/**********************************
* Optical components within cave *
@@ -343,60 +299,14 @@ COMPONENT ac_slit = Slit(
/***************
* Sample area *
***************/
// Beam monitoring before sample
COMPONENT PP_small_sample_focus = PSD_monitor(
filename = "PP_small_sample_focus",
nx = 200, ny = 200,
xwidth = 0.02, yheight = 0.05,
restore_neutron = 1)
WHEN PP_small ==1
AT (0, 0, -0.001) RELATIVE arm_sample
COMPONENT DL_sample_focus = DivLambda_monitor(
filename = "DL_sample_focus",
nL = 200, nh = 200,
xwidth = 0.02, yheight = 0.06,
maxdiv_h = 3,
Lmin = lambda_start, Lmax = lambda_end,
restore_neutron = 1)
AT (0, 0, -0.001) RELATIVE arm_sample
/* monitor at sample position for footprint */
COMPONENT PL_sample = Monitor_nD(
filename = "PL_sample",
options = "lambda limits=[0 35] bins=350 x limits=[-0.02 0.02] bins 160",
xwidth = 0.04, yheight = 0.02, restore_neutron=1)
WHEN sample!=4
AT (0, 0, -0.001) RELATIVE arm_sample
ROTATED (0, -90, 0) RELATIVE arm_sample
COMPONENT PSD_sample = PSD_monitor(
filename = "PSD_sample",
nx = 100, ny = 100,
xwidth = 4*sample_length, yheight = 4*sample_height,
restore_neutron = 1)
WHEN sample!=4
AT (0, 0, -0.001) RELATIVE arm_sample
ROTATED (0, -90, 0) RELATIVE arm_sample
COMPONENT tof_sample = Monitor_nD(
filename = "tof_sample",
options = "xdiv limits=[-0.75 0.75] bins=15 ydiv limits=[-2.0 2.0] bins=40 time limits=[0 0.6] bins=600 lambda limits=[0 35] bins=350",
xwidth=sample_length, yheight = sample_height)
filename = "tof_sample", user1=p_int,
options = "x limits=[-0.025 0.025] bins=1000 y limits=[-0.025 0.025] bins=1000 xdiv limits=[-0.75 0.75] bins=150 ydiv limits=[-2.0 2.0] bins=400 time limits=[0 0.6] bins=6000 lambda limits=[0 35] bins=3500 sy limits=[-1 1] bins=2000 user1 list all",
xwidth=0.05, yheight = 0.05)
WHEN sample==4
AT (0, 0, -0.0005) RELATIVE arm_sample_beam
ROTATED (0, 0, 0) RELATIVE arm_sample_beam
/* monitor for full beam at sample position ( n/cm²/s ) */
COMPONENT PP_small_sample = PSD_monitor(
filename = "PSD_sample_perp",
nx = 200, ny = 200,
xwidth = 0.05, yheight = 0.05,
restore_neutron = 1)
WHEN sample==4
AT (0, 0, 0) RELATIVE arm_sample_beam
ROTATED (0, 0, 0) RELATIVE arm_sample_beam
/* NiTi multilayer sample */
COMPONENT sample = Mirror(
@@ -439,55 +349,67 @@ COMPONENT si_sample = Mirror(
ROTATED (0, 90, 0) RELATIVE arm_sample
COMPONENT arm_analyzer = Arm()
/* Rotate spin 90 degrees*/
COMPONENT arm_hack = Arm()
AT (0, 0, 0) RELATIVE arm_sample
EXTEND %{
double tmp = sy;
sy = sx;
sx = -tmp;
%}
/* flipped arm detector to position analyzer correctly for beam path 1*/
COMPONENT arm_detector2 = Arm()
AT (0, 0, 0) RELATIVE arm_detector
ROTATED (-selene_theta+(Theta1_analyzer1-Theta2_analyzer1)/2.0, 0, 0) RELATIVE arm_detector
ROTATED (0,0,180) RELATIVE arm_detector
COMPONENT arm_analyzer = Arm()
AT (0, 0, 0) RELATIVE arm_detector2
ROTATED (-selene_theta+(Theta1_analyzer1-Theta2_analyzer1)/2.0, 0, 0) RELATIVE arm_detector2
COMPONENT virtual_analyzer_flipper = Arm() // Gone -> Pol_SF_ideal(ny=1, xwidth=1, yheight=1, zdepth=0.001)
WHEN (enable_analyzer>2)
AT (0, 0, analyzer_flipper_start) RELATIVE arm_analyzer
ROTATED (0,0,0.0) RELATIVE arm_analyzer
EXTEND %{
if(INSTRUMENT_GETPAR(enable_analyzer)>2) sx *=-1;
%}
COMPONENT arm_analyzer2 = Arm()
AT (0, 0, 0) RELATIVE arm_detector
ROTATED (-selene_theta+(Theta1_analyzer2-Theta2_analyzer2)/2.0, 0, 0) RELATIVE arm_detector
AT (0, 0, 0) RELATIVE arm_detector2
ROTATED (-selene_theta+(Theta1_analyzer2-Theta2_analyzer2)/2.0, 0, 0) RELATIVE arm_detector2
/* polarization analyser */
COMPONENT analyzer1 = Polariser(nIncRefr=1, d_substrate = 5e-4, reflect_d=0, reflect_u=0, lin=analyzer1_start, length=analyzer1_length,
COMPONENT analyzer1 = Polariser(enable_ref=1, d_substrate = 5e-4, reflect_d=0, reflect_u=0, lin=analyzer1_start, length=analyzer1_length,
delta_theta=(Theta1_analyzer1+Theta2_analyzer1+0.05)*PI/180.0, h2=0.14, h1=0.05, abs_ref=0,
m_u=5.5, m_d=0.45, both_coated=0, alpha=2.3, W = 0.0014)
WHEN enable_analyzer
AT (0, 0.0, analyzer1_start) RELATIVE arm_analyzer
ROTATED (0,0,0.0) RELATIVE arm_analyzer
COMPONENT analyzer2 = Polariser(nIncRefr=1, d_substrate = 5e-4, reflect_d=0, reflect_u=0, lin=analyzer2_start, length=analyzer2_length,
COMPONENT analyzer2 = Polariser(enable_ref=1, d_substrate = 5e-4, reflect_d=0, reflect_u=0, lin=analyzer2_start, length=analyzer2_length,
delta_theta=(Theta1_analyzer2+Theta2_analyzer2+0.05)*PI/180.0, h2=0.2, h1=0.05, abs_out=0,
m_u=5.0, m_d=0.65, both_coated=1, alpha=2.3, W = 0.0014)
WHEN enable_analyzer==2
AT (0, 0.0, analyzer2_start) RELATIVE arm_analyzer
ROTATED (0,0,0.0) RELATIVE arm_analyzer
AT (0, 0.0, analyzer2_start) RELATIVE arm_analyzer2
ROTATED (0,0,0.0) RELATIVE arm_analyzer2
/* detector */
// Beam monitoring at detector position
COMPONENT PP_large_detector = PSD_monitor(
filename = "PP_large_detector",
nx = 200, ny = 200,
xmin = -1.8, xmax = 1.8, ymin = -1.8, ymax = 1.8,
restore_neutron = 1)
AT (0, 0, total_length+detector_arm-0.001) RELATIVE ISCS
// Detector with 0.5x1mm² resolution
COMPONENT PP_detector = PSD_monitor(
filename = "PP_detector",
nx = 1000, ny = 500,
xwidth = 0.5, yheight = 0.5,
restore_neutron = 1)
AT (0, 0, detector_arm) RELATIVE arm_detector
COMPONENT tof_detector = Monitor_nD(
filename = "tof_detector",
options = "x limits=[-0.25 0.25] bins=1000 y limits=[-0.5 0.5] bins=1000 time limits=[0 0.6] bins=6000 lambda limits=[0 35] bins=3500 sx limits=[-1 1] bins=1000 sy limits=[-1 1] bins=1000, list all",
filename = "tof_detector", user1=p_int,
options = "x limits=[-0.25 0.25] bins=1000 y limits=[-0.5 0.5] bins=1000 time limits=[0 0.6] bins=6000 lambda limits=[0 35] bins=3500 sx sy sz user1, list all",
xwidth = 0.5, yheight = 1.0)
WHEN sample!=4
AT (0, 0, detector_arm+0.00001) RELATIVE arm_detector
/*COMPONENT detector = PSD_TOF_monitor(xwidth=0.5, yheight=0.5,
nx=1000, ny=250, tmin=35000, tmax=115000, nt=160, filename="detector")
AT (0, 0, detector_arm) RELATIVE arm_detector
ROTATED (0, 0, 0) RELATIVE arm_detector
/***********************************************************************/
+17 -17
View File
@@ -144,7 +144,7 @@ COMPONENT NBOA_window=Al_window(thickness=NBOA_Al_entrance_length)
AT (0,0,NBOA_Al_entrance_start) RELATIVE ISCS
COMPONENT NBOA_side = AbsorberFixed(xmin=NBOA_side_x, xmax=NBOA_side_x+0.012,
COMPONENT NBOA_side = Absorber(xmin=NBOA_side_x, xmax=NBOA_side_x+0.012,
ymin=-0.2, ymax=0.2,
zmin=0.0, zmax=E02_01_01_Cu_length+E02_01_01_length+E02_01_02_length+E02_01_03_length)
AT (0, 0, E02_01_01_Cu_start) RELATIVE ISCS
@@ -181,7 +181,7 @@ COMPONENT NBOA_Cu_collimator=Guide_gravity(
* feeder neutron guide *
************************/
COMPONENT E02_01_011_wedge = AbsorberFixed(xmin=-0.1, xmax=0.0,
COMPONENT E02_01_011_wedge = Absorber(xmin=-0.1, xmax=0.0,
ymin=-E02_01_011_wh/2.0, ymax=E02_01_011_wh/2.0,
zmin=0.0, zmax=E02_01_01_length-0.5)
AT (0, 0, E02_01_01_start) RELATIVE ISCS
@@ -190,7 +190,7 @@ COMPONENT E02_01_011_wedge = AbsorberFixed(xmin=-0.1, xmax=0.0,
ALLOW_BACKPROP;
PROP_Z0;
%}
COMPONENT E02_01_012_wedge = AbsorberFixed(xmin=-0.1, xmax=0.0,
COMPONENT E02_01_012_wedge = Absorber(xmin=-0.1, xmax=0.0,
ymin=-E02_01_012_wh/2.0, ymax=E02_01_012_wh/2.0,
zmin=0.5-E02_01_01_Cu_length, zmax=E02_01_01_length)
AT (0, 0, E02_01_01_start) RELATIVE ISCS
@@ -210,7 +210,7 @@ COMPONENT E02_01_01 = Elliptic_guide_gravity(
enableGravity=enable_gravity)
AT (0, 0, E02_01_01_start) RELATIVE ISCS
COMPONENT E02_01_021_wedge = AbsorberFixed(xmin=-0.1, xmax=0.0,
COMPONENT E02_01_021_wedge = Absorber(xmin=-0.1, xmax=0.0,
ymin=-E02_01_021_wh/2.0, ymax=E02_01_021_wh/2.0,
zmin=0.0, zmax=E02_01_02_length/2.0)
AT (0, 0, E02_01_02_start) RELATIVE ISCS
@@ -219,7 +219,7 @@ COMPONENT E02_01_021_wedge = AbsorberFixed(xmin=-0.1, xmax=0.0,
ALLOW_BACKPROP;
PROP_Z0;
%}
COMPONENT E02_01_022_wedge = AbsorberFixed(xmin=-0.1, xmax=0.0,
COMPONENT E02_01_022_wedge = Absorber(xmin=-0.1, xmax=0.0,
ymin=-E02_01_022_wh/2.0, ymax=E02_01_022_wh/2.0,
zmin=E02_01_02_length/2.0, zmax=E02_01_02_length)
AT (0, 0, E02_01_02_start) RELATIVE ISCS
@@ -238,7 +238,7 @@ COMPONENT E02_01_02 = Elliptic_guide_gravity(
enableGravity=enable_gravity)
AT (0, 0, E02_01_02_start) RELATIVE ISCS
COMPONENT E02_01_031_wedge = AbsorberFixed(xmin=-0.1, xmax=0.0,
COMPONENT E02_01_031_wedge = Absorber(xmin=-0.1, xmax=0.0,
ymin=-E02_01_031_wh/2.0, ymax=E02_01_031_wh/2.0,
zmin=0.0, zmax=0.5)
AT (0, 0, E02_01_03_start) RELATIVE ISCS
@@ -247,7 +247,7 @@ COMPONENT E02_01_031_wedge = AbsorberFixed(xmin=-0.1, xmax=0.0,
ALLOW_BACKPROP;
PROP_Z0;
%}
COMPONENT E02_01_032_wedge = AbsorberFixed(xmin=-0.1, xmax=0.0,
COMPONENT E02_01_032_wedge = Absorber(xmin=-0.1, xmax=0.0,
ymin=-E02_01_032_wh/2.0, ymax=E02_01_032_wh/2.0,
zmin=0.5, zmax=1.0)
AT (0, 0, E02_01_03_start) RELATIVE ISCS
@@ -256,7 +256,7 @@ COMPONENT E02_01_032_wedge = AbsorberFixed(xmin=-0.1, xmax=0.0,
ALLOW_BACKPROP;
PROP_Z0;
%}
COMPONENT E02_01_033_wedge = AbsorberFixed(xmin=-0.1, xmax=0.0,
COMPONENT E02_01_033_wedge = Absorber(xmin=-0.1, xmax=0.0,
ymin=-E02_01_033_wh/2.0, ymax=E02_01_033_wh/2.0,
zmin=1.0, zmax=E02_01_03_length)
AT (0, 0, E02_01_03_start) RELATIVE ISCS
@@ -299,7 +299,7 @@ COMPONENT NFGA_Al_window_in=Al_window(
WHEN enable_windows
AT (0,0,NFGA_Al_start) RELATIVE ISCS
COMPONENT NFGA_side1 = AbsorberFixed(xmin=NFGA_side_x1, xmax=NFGA_side_x1+0.008,
COMPONENT NFGA_side1 = Absorber(xmin=NFGA_side_x1, xmax=NFGA_side_x1+0.008,
ymin=-0.2, ymax=0.2,
zmin=0.002, zmax=0.502)
AT (0, 0, E02_02_01_start-0.002) RELATIVE ISCS
@@ -309,7 +309,7 @@ COMPONENT NFGA_side1 = AbsorberFixed(xmin=NFGA_side_x1, xmax=NFGA_side_x1+0.008,
ALLOW_BACKPROP;
PROP_Z0;
%}
COMPONENT NFGA_side2 = AbsorberFixed(xmin=NFGA_side_x2, xmax=NFGA_side_x2+0.008,
COMPONENT NFGA_side2 = Absorber(xmin=NFGA_side_x2, xmax=NFGA_side_x2+0.008,
ymin=-0.2, ymax=0.2,
zmin=0.502, zmax=E02_02_01_length+E02_02_02_length+E02_02_03_length+0.002)
AT (0, 0, E02_02_01_start-0.002) RELATIVE ISCS
@@ -321,7 +321,7 @@ COMPONENT NFGA_side2 = AbsorberFixed(xmin=NFGA_side_x2, xmax=NFGA_side_x2+0.008,
%}
// in-bunker feeder segments
COMPONENT E02_02_011_wedge = AbsorberFixed(xmin=-0.1, xmax=0.0,
COMPONENT E02_02_011_wedge = Absorber(xmin=-0.1, xmax=0.0,
ymin=-E02_02_011_wh/2.0, ymax=E02_02_011_wh/2.0,
zmin=0.0, zmax=E02_02_01_length/2.0)
AT (0, 0, E02_02_01_start) RELATIVE ISCS
@@ -330,7 +330,7 @@ COMPONENT E02_02_011_wedge = AbsorberFixed(xmin=-0.1, xmax=0.0,
ALLOW_BACKPROP;
PROP_Z0;
%}
COMPONENT E02_02_012_wedge = AbsorberFixed(xmin=-0.1, xmax=0.0,
COMPONENT E02_02_012_wedge = Absorber(xmin=-0.1, xmax=0.0,
ymin=-E02_02_012_wh/2.0, ymax=E02_02_012_wh/2.0,
zmin=E02_02_01_length/2.0, zmax=E02_02_01_length)
AT (0, 0, E02_02_01_start) RELATIVE ISCS
@@ -349,7 +349,7 @@ COMPONENT E02_02_01 = Elliptic_guide_gravity(
enableGravity=enable_gravity)
AT (0, 0, E02_02_01_start) RELATIVE ISCS
COMPONENT E02_02_021_wedge = AbsorberFixed(xmin=-0.1, xmax=0.0,
COMPONENT E02_02_021_wedge = Absorber(xmin=-0.1, xmax=0.0,
ymin=-E02_02_021_wh/2.0, ymax=E02_02_021_wh/2.0,
zmin=0.0, zmax=E02_02_02_length/2.0)
AT (0, 0, E02_02_02_start) RELATIVE ISCS
@@ -358,7 +358,7 @@ COMPONENT E02_02_021_wedge = AbsorberFixed(xmin=-0.1, xmax=0.0,
ALLOW_BACKPROP;
PROP_Z0;
%}
COMPONENT E02_02_022_wedge = AbsorberFixed(xmin=-0.1, xmax=0.0,
COMPONENT E02_02_022_wedge = Absorber(xmin=-0.1, xmax=0.0,
ymin=-E02_02_022_wh/2.0, ymax=E02_02_022_wh/2.0,
zmin=E02_02_02_length/2.0, zmax=E02_02_02_length)
AT (0, 0, E02_02_02_start) RELATIVE ISCS
@@ -377,7 +377,7 @@ COMPONENT E02_02_02 = Elliptic_guide_gravity(
enableGravity=enable_gravity)
AT (0, 0, E02_02_02_start) RELATIVE ISCS
COMPONENT E02_02_031_wedge = AbsorberFixed(xmin=-0.1, xmax=0.0,
COMPONENT E02_02_031_wedge = Absorber(xmin=-0.1, xmax=0.0,
ymin=-E02_02_031_wh/2.0, ymax=E02_02_031_wh/2.0,
zmin=0.0, zmax=0.5)
AT (0, 0, E02_02_03_start) RELATIVE ISCS
@@ -386,7 +386,7 @@ COMPONENT E02_02_031_wedge = AbsorberFixed(xmin=-0.1, xmax=0.0,
ALLOW_BACKPROP;
PROP_Z0;
%}
COMPONENT E02_02_032_wedge = AbsorberFixed(xmin=-0.1, xmax=0.0,
COMPONENT E02_02_032_wedge = Absorber(xmin=-0.1, xmax=0.0,
ymin=-E02_02_032_wh/2.0, ymax=E02_02_032_wh/2.0,
zmin=0.5, zmax=1.0)
AT (0, 0, E02_02_03_start) RELATIVE ISCS
@@ -395,7 +395,7 @@ COMPONENT E02_02_032_wedge = AbsorberFixed(xmin=-0.1, xmax=0.0,
ALLOW_BACKPROP;
PROP_Z0;
%}
COMPONENT E02_02_033_wedge = AbsorberFixed(xmin=-0.1, xmax=0.0,
COMPONENT E02_02_033_wedge = Absorber(xmin=-0.1, xmax=0.0,
ymin=-E02_02_033_wh/2.0, ymax=E02_02_033_wh/2.0,
zmin=1.0, zmax=E02_02_03_length)
AT (0, 0, E02_02_03_start) RELATIVE ISCS
+244
View File
@@ -0,0 +1,244 @@
/*******************************************************************************
* McStas instrument definition URL=http://www.mcstas.org
*
* Instrument: Estia_mf
*
* %Identification
* Written by: Artur Glavic (artur.glavic@psi.ch); Jochen Stahn (jochen.stahn@psi.ch); Christine Klauser (christine.klauser@psi.ch)
* Date: 01. 03. 2018
* Origin: PSI
* Release: McStas 2.4.1
* Version: 1.0
* %INSTRUMENT_SITE: ESS (E02)
*
* Estia is a vertical sample, focusing reflectometer for small sample
*
* %Description
* These are the components near the middle focus between the Selene guide 1 and 2.
* This file is not intended for stand alone use but is included in the Estia instrument model.
* The removable components, however, allow for a stand alone use for debugging purpose.
*
* %Parameters
* enable_gravity: [0/1] Use gravity in elliptical guide model.
*
* %End
*******************************************************************************/
DEFINE INSTRUMENT Estia_selene(int enable_gravity=0)
DECLARE
%{
// Polarizer parameters used for the detailed design
double pol1_xstart = 1.007; // [m] up-stream of focus point
double pol1_ystart = 0.0064; // [m] (relative to c-axis)
double pol1_xend = 0.365; // [m] up-stream of focus point
double pol1_yend = 0.013; // [m]
double pol1_gamma = 1.65; // [deg] optimal incident angle for beams passing through focus point
double pol1_m = 5.08;
double pol1_lin_xstart = 1.007; //[m]
double pol1_lin_ystart = 0.0064; //[m]
double pol1_lin_xend = 0.812; //[m]
double pol1_lin_yend = -0.0014; //[m]
double pol1_lin_m = 5.5;
double pol2_xstart = 1.286; // [m] up-stream of focus point
double pol2_ystart = 0.009; // [m] (relative to c-axis)
double pol2_xend = 0.500; // [m] up-stream of focus point
double pol2_yend = 0.0175; // [m]
double pol2_gamma = 1.70; // [deg] optimal incident angle for beams passing through focus point
double pol2_m = 5.08;
double pol2_lin_xstart = 1.286; //[m]
double pol2_lin_ystart = 0.009; //[m]
double pol2_lin_xend = 1.093; //[m]
double pol2_lin_yend = 0.0014; //[m]
double pol2_lin_m = 5.08;
double pol_hclose = 0.12; // [m] height of entrance
double pol_hfar = 0.12; // [m] height of exit
// parameters for polarizer component calculated from the above
double pol1_length, pol2_length, pol1_lin_length, pol2_lin_length;
double Theta_pol1, Theta_pol2, Theta_rot1, Theta_rot2;
double Theta1_pol1, Theta1_pol2, Theta2_pol1, Theta2_pol2;
double pol1_lin_rot, pol2_lin_rot;
%}
INITIALIZE
%{
pol1_length = pol1_xstart-pol1_xend;
pol2_length = pol2_xstart-pol2_xend;
pol1_lin_length = pol1_lin_xstart-pol1_lin_xend;
pol2_lin_length = pol2_lin_xstart-pol2_lin_xend;
Theta1_pol1=atan2(pol1_ystart, pol1_xstart); // [rad] angle between c-axis and entrance point
Theta2_pol1=atan2(pol1_yend, pol1_xend); // [rad] angle between c-axis and exit point
Theta1_pol2=atan2(pol2_ystart, pol2_xstart); // [rad] angle between c-axis and entrance point
Theta2_pol2=atan2(pol2_yend, pol2_xend); // [rad] angle between c-axis and exit point
Theta_pol1=(Theta2_pol1-Theta1_pol1)*180.0/PI; // [deg] full covered divergence angle
Theta_pol2=(Theta2_pol2-Theta1_pol2)*180.0/PI; // [deg]
Theta_rot1=(Theta1_pol1+Theta2_pol1)/2.0*180.0/PI; // [deg] rotation of center of polarizer
Theta_rot2=(Theta1_pol2+Theta2_pol2)/2.0*180.0/PI; // [deg]
pol1_lin_rot=atan2(pol1_lin_ystart-pol1_lin_yend,
pol1_lin_xstart-pol1_lin_xend)*180.0/PI; // [deg] rotation angle of polarizer
pol2_lin_rot=atan2(pol2_lin_ystart-pol2_lin_yend,
pol2_lin_xstart-pol2_lin_xend)*180.0/PI; // [deg] rotation angle of polarizer
printf(" Polarizer 1 angle = %.2f deg\n", Theta_rot1);
printf(" Polarizer 1 divergence = %.2f deg\n", Theta_pol1);
printf(" Polarizer 1 distance = %.1f mm\n", pol1_xstart*1000.0);
printf(" Polarizer 1 length = %.1f mm\n", pol1_length*1000.0);
printf(" PolStraight 1 angle = %.2f deg\n", pol1_lin_rot);
printf(" Polarizer 2 angle = %.2f deg\n", Theta_rot1);
printf(" Polarizer 2 divergence = %.2f deg\n", Theta_pol2);
printf(" Polarizer 2 distance = %.1f mm\n", pol2_xstart*1000.0);
printf(" Polarizer 2 length = %.1f mm\n", pol2_length*1000.0);
printf(" PolStraight 2 angle = %.2f deg\n", pol2_lin_rot);
%}
TRACE
REMOVABLE COMPONENT origin = Progress_bar()
AT (0,0,0) ABSOLUTE
REMOVABLE COMPONENT ISCS = Arm()
AT (0, 0, 0) RELATIVE origin
REMOVABLE COMPONENT arm_selene1 = Arm()
AT (0, 0, 0) RELATIVE ISCS
REMOVABLE COMPONENT arm_selene2 = Arm()
AT (0, 0, 2*selene_c) RELATIVE ISCS
REMOVABLE COMPONENT source = Moderator(radius = 0.001, focus_xw = selene_entry+0.002, focus_yh = selene_entry+0.002,
target_index=3, Emin=1, Emax=15)
AT (0,0,0) RELATIVE ISCS
ROTATED (-1.25, 1.25, 0) RELATIVE ISCS
/**************************************
* Middle focus between Selene guides *
**************************************/
COMPONENT arm_polref = Arm()
AT (0, 0, 2*selene_c) RELATIVE arm_selene1
ROTATED (selene_theta, -selene_theta, 0) RELATIVE ISCS
COMPONENT arm_polarizer = Arm()
AT (0, 0, 2*selene_c) RELATIVE arm_selene1
ROTATED (0, 0, 0) RELATIVE arm_selene1
COMPONENT arm_pol1 = Arm()
AT (0, 0, 0) RELATIVE arm_polarizer
ROTATED (0, -Theta_rot1, 0) RELATIVE arm_polarizer
COMPONENT arm_pol2 = Arm()
AT (0, 0, 0) RELATIVE arm_polarizer
ROTATED (0, -Theta_rot2, 0) RELATIVE arm_polarizer
COMPONENT polarizer2_lin = Pol_mirror(zwidth = pol2_lin_length, yheight = pol_hfar, p_reflect=0,
rUpFunc=TableReflecFunc, rUpPar="SNAG_m5p5_Tup.ref",
rDownFunc=TableReflecFunc, rDownPar="SNAG_m5p5_Tdown.ref", useTables=1)
WHEN (enable_polarizer>0)
AT ((pol2_lin_ystart+pol2_lin_yend)/2.0, 0, -(pol2_lin_xstart+pol2_lin_xend)/2.0) RELATIVE arm_polarizer
ROTATED (0,-pol2_lin_rot,0.0) RELATIVE arm_polarizer
GROUP polarizer2_set
EXTEND
%{
if (SCATTERED==2) {
p_int +=2;
}
%}
COMPONENT polarizer2 = Polariser(lin=-pol2_xstart, length=pol2_length,
enable_ref=1, abs_ref=1, abs_out=1, both_coated=1,
d_substrate = 5e-4, T_loss=4.0e3,
m_u=pol2_m, m_d=0.6, m_residual=0.55,
alpha=2.3, W = 0.0014, reflect_d=0, reflect_u=0,
delta_theta=Theta_pol2*PI/180.0,
h2=pol_hfar, h1=pol_hclose)
WHEN (enable_polarizer>0)
AT (0, 0, -pol2_xstart) RELATIVE arm_pol2
ROTATED (0,0,90.0) RELATIVE arm_pol2
GROUP polarizer2_set
EXTEND
%{
if (SCATTERED) {
ALLOW_BACKPROP;
PROP_Z0;
p_int +=1;
} else ABSORB;
%}
COMPONENT replacement_pol2_set = Slit(xwidth=0.2, yheight=0.2)
WHEN (enable_polarizer==0)
AT (0, 0, -0.2) RELATIVE arm_polarizer
GROUP polarizer2_set
COMPONENT polarizer1_lin = Pol_mirror(zwidth = pol1_lin_length, yheight = pol_hfar, p_reflect=0,
rUpFunc=TableReflecFunc, rUpPar="SNAG_m5p5_Tup.ref",
rDownFunc=TableReflecFunc, rDownPar="SNAG_m5p5_Tdown.ref", useTables=1)
WHEN (enable_polarizer>0)
AT ((pol1_lin_ystart+pol1_lin_yend)/2.0, 0, -(pol1_lin_xstart+pol1_lin_xend)/2.0) RELATIVE arm_polarizer
ROTATED (0,-pol1_lin_rot,0.0) RELATIVE arm_polarizer
GROUP polarizer1_set
EXTEND
%{
if (SCATTERED==2) {
p_int +=8;
}
%}
COMPONENT polarizer1 = Polariser(lin=-pol1_xstart, length=pol1_length,
enable_ref=1, abs_ref=1, abs_out=1, both_coated=1,
d_substrate = 5e-4, T_loss=4.0e3,
m_u=pol1_m, m_d=0.6, m_residual=0.55,
alpha=2.3, W = 0.0014, reflect_d=0, reflect_u=0,
delta_theta=Theta_pol1*PI/180.0,
h2=pol_hfar, h1=pol_hclose)
WHEN (enable_polarizer>0)
AT (0, 0, -pol1_xstart) RELATIVE arm_pol1
ROTATED (0,0,90.0) RELATIVE arm_pol1
GROUP polarizer1_set
EXTEND
%{
if (SCATTERED) {
p_int +=4;
}
%}
COMPONENT replacement_pol1_set = Slit(xwidth=0.2, yheight=0.2)
WHEN (enable_polarizer==0)
AT (0, 0, -0.1) RELATIVE arm_polarizer
GROUP polarizer1_set
COMPONENT mf_divergence_h = DivPos_monitor(nb=21, ndiv=41, filename="mf_hordiv",
xmin=-0.02, xmax=0.02, ymin=-0.02, ymax=0.02, maxdiv=2.0)
AT (0, 0, 0) RELATIVE arm_polarizer
COMPONENT mf_divergence_v = DivPos_monitor(nb=21, ndiv=41, filename="mf_verdiv",
xmin=-0.02, xmax=0.02, ymin=-0.02, ymax=0.02, maxdiv=2.0, vertical=1)
AT (0, 0, 0) RELATIVE arm_polarizer
COMPONENT virtual_flipper = Pol_SF_ideal(ny=1, xwidth=0.1, yheight=0.1, zdepth=0.001)
WHEN (enable_polarizer>2)
AT (0, 0, 0) RELATIVE arm_polarizer
FINALLY
%{
%}
END
+204
View File
@@ -0,0 +1,204 @@
/*******************************************************************************
* McStas instrument definition URL=http://www.mcstas.org
*
* Instrument: ESS_reflectometer_Estia
*
* %Identification
* Written by: Artur Glavic (artur.glavic@psi.ch); Jochen Stahn (jochen.stahn@psi.ch); Christine Klauser (christine.klauser@psi.ch)
* Date: 01. 03. 2018
* Origin: PSI
* Release: McStas 2.4.1
* Version: 1.0
* %INSTRUMENT_SITE: ESS (E02)
*
* Estia is a vertical sample, focusing reflectometer for small sample
*
* %Description
* The instrument consists of a two part elliptical feeder that focuses onto
* a slit mask that defines the sample footprint, called virtual source.
* The virtual source is imaged on the sample position with a Selene
* type neutron guide (two ellipses).
* This version of the instrument, used in the ESS butterfly moderator,
* has two vertical Selene guide systems that share the same focal points.
* The two guides are only implemented in the feeder as the selene mirrors
* are part of the upgrade program.
*
* %Parameters
* omegaa: [deg] sample rotation omega
* sample: flag to switch between (0) sample Ni/Ti-multilayer
* (1) reference ( R(q_z) = 1 for all q_z )
* (2) Ni-film on silicon
* (3) Natural SiO2 on silicon
* (4) monitor instead of sample (normal to beam)
* sample_length: [m] Size of sample in beam direction, also controls virtual source opening
* sample_height: [m] Size of sample in vertical direction, also controls virtual source opening
*
* operationmode: operation mode (0) high-intensity specular reflectivity
* (1) almost conventional TOF
* over_illumination: [m] Extra opening of virtual source compared to sample size
* theta_resolution: Delta theta / theta adjusted with the slit (operationmode=1)
*
* lambda_start: [A] Beginning of simulated wavelength range
* lambda_end: [A] End of simulated wavelength range
* enable_gravity: [0/1] Use gravity in elliptical guide model.
* enable_chopper: [0/1/2/3] Activate chopper component if !=0. Numbers larger than
* 1 define the number of pulses per chopper opening.
* (pulses skipped=enable_chopper-1)
*
* %End
*******************************************************************************/
DEFINE INSTRUMENT ESS_reflectometer_Estia
(double lambda_start = 3.0, double lambda_end = 12.0,
int enable_gravity=1, int enable_windows=1, direct_beam=0,
double source_power = 5, foil_thickness=0.00001)
DECLARE
%{
/* Geometrical parameters from CAD model of Estia (ESS-0050413)
* TCS coordinate and directional rotation first focus point
* refered to as focus_moderator_y_rot
*/
//TCS position of ISCS: (110,-105,137)
double iscs_x=0.0199717;
double iscs_y=0.0;
double iscs_z=0.0547095; // ISCS in McStas coordinates using mod-view-opt.instr
double iscs_rot_y=0.4; // ISCS is at TCS-35.6 degree, McStas at -36 degree
double iscs_rot_x=0.7; // downward tilt of Estia axis
// Selene 1 geometry parameters (optics parameters in Estia_selene1.instr)
double selene1_foot1 = 4.20; // distance of first foot to VS focus
double selene1_foot2 = 7.00; // distance of second foot to VS focus
double selene1_center;
double selene1_shift;
double selene1_rot;
// Selene 2
/* general instrument geometry parameters */
double total_length = 35.0 ; // m distance moderator-sample (2*c_feeder+4*c_Selene)
double detector_arm = 4.0 ; // m distance sample-detector
/* chopper parameters */
double chopper_diameter = 0.7;
double chopper_pos = 10.895; // m distance source-chopper
double chopper_phase ; // deg phase between pulse and chopper opening
double chopper_open ; // deg of opening angle in the chopper
double pulse_zero = 0.00175; // ms intensity weighted average time of emitted neutron pulse
double opening_time = 0.0015; // ms time to reach full intensity on detector, used to adjust phase to get full intensity at beginning of selected band
double chopper_freq = 14.0 ; // Hz chopper frequency
double slit_distance = 1.775; //m, distance slit to sample
/* derived quantities */
double selene_theta = 1.25;
%}
INITIALIZE
%{
%}
TRACE
/********************************************************************************
* Initial geometry coordinate axes for important components of the simulation. *
********************************************************************************/
COMPONENT origin = Progress_bar()
AT (0,0,0) ABSOLUTE
/* The ISCS is the instrument coordinate system with directions where y points upwards and z lies on the instrument axes. */
COMPONENT ISCS_rot1 = Arm() // position correctly and rotate around z-axis (x points in beam direction)
AT (iscs_x,iscs_y,iscs_z) RELATIVE origin
ROTATED (0,iscs_rot_y,0) RELATIVE origin
COMPONENT ISCS = Arm() // rotate around y-axis (slight downward tilt)
AT (0,0,0) RELATIVE ISCS_rot1
ROTATED (iscs_rot_x,0,0) RELATIVE ISCS_rot1
/***********************
* Instrument Skeleton *
***********************/
// Axes system parallel to the c-axis of the feeder
COMPONENT arm_feeder = Arm()
AT (0, 0, 0) RELATIVE ISCS
ROTATED (0, 0, 0) RELATIVE ISCS
COMPONENT arm_selene1_center = Arm()
AT (selene1_shift, 0, 2*NBOA_c+selene1_center) RELATIVE ISCS
ROTATED (0, selene1_rot, 0) RELATIVE ISCS
// Axes starting at focus of feeder parallel to c-axis of Selene guides
COMPONENT arm_selene1 = Arm()
AT (0, 0, -selene1_center) RELATIVE arm_selene1_center
// Axes of the center of usable beam passing through the virtual source
COMPONENT arm_virtual_source_beam = Arm()
AT (0, 0, 0) RELATIVE arm_selene1
ROTATED (0, selene_theta, 0) RELATIVE ISCS
// 90 degree to main beam is the monitors detector position
COMPONENT arm_monitor = Arm()
AT (0, 0, -0.320) RELATIVE arm_virtual_source_beam
ROTATED (-90, 0, 0) RELATIVE arm_virtual_source_beam
/**********
* Source *
**********/
COMPONENT moderator = ESS_butterfly(
sector = "E", beamline = 2, yheight = 0.03, cold_frac = 0.9,
focus_xw = E02_01_01_Cu_in_xmax-E02_01_01_Cu_in_xmin+0.002,
focus_yh = E02_01_01_Cu_in_yheight+0.002,
target_index=4, Lmin = lambda_start, Lmax = lambda_end,
n_pulses = 1, acc_power=source_power)
AT (0, 0, 0) RELATIVE origin
/***************************************
* Geometry of neutron feeder separate *
***************************************/
%include "Estia_feeder.instr"
/****************************************************
* Beam manipulation area around the virtual source *
****************************************************/
/* Absorber to reduce beam to needed size an for shielding purposes (CPC1 in CAD model) */
COMPONENT CPC1_in = Slit(xwidth=0.0303, yheight=0.0792)
AT (0, 0, -0.890) RELATIVE arm_virtual_source_beam
COMPONENT CPC1_monitor = Slit(xwidth=0.016344, yheight=0.05547)
AT (0, 0, -0.3573) RELATIVE arm_virtual_source_beam
// As TOF detector is rectangular, use focus_r size to limit to actual
// detector size of 0.5'' diameter cylinder
COMPONENT Vanadium_Foil = Incoherent(focus_r=0.00635, p_interact=0.9,
xwidth=0.018, yheight=0.055, zdepth=foil_thickness,
sigma_abs=5.08, sigma_inc=5.08, Vc=13.827,
target_index=1)
WHEN direct_beam<2
AT (0, 0, -0.320) RELATIVE arm_virtual_source_beam
ROTATED (34.3, 0, 0) RELATIVE arm_virtual_source_beam
COMPONENT tof_monitor = TOFLambda_monitor(
filename = "monitor",
tmin=0, tmax=120000, nt=1200,
Lmin=0,Lmax=35,nL=350,
xwidth = 0.025, yheight = 0.025)
WHEN direct_beam==0
AT (0, 0, 0.1) RELATIVE arm_monitor
COMPONENT tof_direct = TOFLambda_monitor(
filename = "VS",
tmin=0, tmax=120000, nt=1200,
Lmin=0,Lmax=35,nL=350,
xwidth = 0.05, yheight = 0.05)
WHEN direct_beam>0
AT (0, 0, 0.08) RELATIVE arm_virtual_source_beam
FINALLY
%{
%}
END
-240
View File
@@ -1,240 +0,0 @@
/*******************************************************************************
* McStas instrument definition URL=http://www.mcstas.org
*
* Instrument: Estia_selene
*
* %Identification
* Written by: Artur Glavic (artur.glavic@psi.ch); Jochen Stahn (jochen.stahn@psi.ch); Christine Klauser (christine.klauser@psi.ch)
* Date: 01. 03. 2018
* Origin: PSI
* Release: McStas 2.4.1
* Version: 1.0
* %INSTRUMENT_SITE: ESS (E02)
*
* Estia is a vertical sample, focusing reflectometer for small sample
*
* %Description
* This is the beam extraction part of the instrument, it is not intended for use but is
* included in the Estia instrument model. The removable components, however,
* allow for a stand alone use for debugging purpose.
*
* %Parameters
* enable_gravity: [0/1] Use gravity in elliptical guide model.
*
* %End
*******************************************************************************/
DEFINE INSTRUMENT Estia_selene(int enable_gravity=0)
DECLARE
%{
// Selene 1 (Will follow later with detailed parameters)
// Selene 2
/* selene parameters */
double selene_xi = 0.6 ; // length of neutron guide over 2*c of ellipses, xi
double selene_b = 0.1047; // m selene short half axis c
double selene_c = 6.0000; // m selene focus half distance c
double selene_coating = 3.8; // m-value of selene coating
/* other variables */
double selene_entry ;
double selene_length ;
double selene_distance ;
double polarizer_max_width = 0.0015; // Maximum sample height to be covered by te polarization analyzers
double polarizer_start = 0.3; // Distance to sample to start the first analyzer
double polarizer_length = 0.6; // length of first analyzer
double Theta1_polarizer, Theta2_polarizer, dist_pol_vfocus; // quantities calculated out of values above and 1.5 degree covered divergence
%}
INITIALIZE
%{
/* derived parameters for the selene guides */
selene_length = 2.0 * selene_c * selene_xi;
selene_distance = ( 1.0 - selene_xi ) * selene_c;
selene_entry = selene_b * sqrt( 1.0 - pow(selene_xi,2.0) );
dist_pol_vfocus=polarizer_max_width/2.0/tan(1.5*PI/360.0); // the virtual focus point infront of the actual sample focus where the beams furthest out meet
Theta1_polarizer=atan((dist_pol_vfocus+polarizer_start)/dist_pol_vfocus*polarizer_max_width/2.0/polarizer_start)*180.0/PI;
Theta2_polarizer=atan((dist_pol_vfocus+polarizer_start+polarizer_length)/dist_pol_vfocus*polarizer_max_width/2.0/(polarizer_start+polarizer_length))*180.0/PI;
%}
TRACE
REMOVABLE COMPONENT origin = Progress_bar()
AT (0,0,0) ABSOLUTE
REMOVABLE COMPONENT ISCS = Arm()
AT (0, 0, 0) RELATIVE origin
REMOVABLE COMPONENT arm_selene1 = Arm()
AT (0, 0, 0) RELATIVE ISCS
REMOVABLE COMPONENT arm_selene2 = Arm()
AT (0, 0, 2*selene_c) RELATIVE ISCS
REMOVABLE COMPONENT source = Moderator(radius = 0.001, focus_xw = selene_entry+0.002, focus_yh = selene_entry+0.002,
target_index=3, Emin=1, Emax=15)
AT (0,0,0) RELATIVE ISCS
ROTATED (-1.25, 1.25, 0) RELATIVE ISCS
/**************************
* Selene 1 neutron guide *
**************************/
/* Absorber to cut direct view beam (Copper in CAD model) */
COMPONENT slit_before_selene_guide_1 = Slit(
xmin = 0, xmax = selene_entry+0.001,
ymin = 0, ymax = selene_entry+0.005)
AT (0, 0, selene_distance-0.01) RELATIVE arm_selene1
COMPONENT block_before_selene_guide_1 = Absorber(
xmin=-1, xmax=1,
ymin=-selene_entry/4.0+0.005,
ymax=selene_entry/4.0-0.005,
zmin=0.0, zmax=0.001)
AT (0, 0, selene_distance-0.0095) RELATIVE arm_selene1
/* Selene 1 elliptic guide */
COMPONENT selene_guide_1 = Elliptic_guide_gravity(
l=selene_length, dimensionsAt = "mid",
linyh = selene_distance, loutyh= selene_distance,
linxw = selene_distance, loutxw= selene_distance,
xwidth=selene_b*2, yheight=selene_b*2,
mright=0, mleft=selene_coating, mtop=selene_coating, mbottom=0,
enableGravity=enable_gravity)
AT (0, 0, selene_distance) RELATIVE arm_selene1
/* Absorber to cut direct view beam (Copper in CAD model) */
COMPONENT slit_after_selene_guide_1 = Slit(
xmin = 0, xmax = selene_entry+0.005,
ymin = 0, ymax = selene_entry+0.01)
AT (0, 0, 2*selene_c-selene_distance+0.001) RELATIVE arm_selene1
COMPONENT block_after_selene_guide_1 = Absorber(
xmin=-1, xmax=1,
ymin=-selene_entry/4.0+0.015,
ymax=selene_entry/4.0-0.015,
zmin=0.0, zmax=0.001)
AT (0, 0, 2*selene_c-selene_distance+0.0015) RELATIVE arm_selene1
/**************************************
* Middle focus between Selene guides *
**************************************/
COMPONENT arm_polarizer = Arm()
AT (0, 0, 2*selene_c) RELATIVE arm_selene1
ROTATED (selene_theta, -selene_theta, 0) RELATIVE ISCS
COMPONENT polarizer1 = Polariser(nIncRefr=1, d_substrate = 5e-4, reflect_d=0, reflect_u=0, lin=-(polarizer_start+polarizer_length), length=polarizer_length,
delta_theta=(Theta1_polarizer+Theta2_polarizer)*PI/180.0, h2=0.1, h1=0.05, abs_ref=1, m_u=4.0, m_d=0.65, both_coated=1, alpha=2.3, W = 0.0014)
WHEN enable_polarizer
AT (0, 0, -(polarizer_start+polarizer_length)) RELATIVE arm_polarizer
ROTATED (0,0,90.0) RELATIVE arm_polarizer
COMPONENT polmon_mfocus=PolLambda_monitor(Lmin=lambda_start, Lmax=lambda_end, nL=30, npol=101, my=1, filename="polmon_mfocus")
WHEN enable_polarizer
AT (0, 0, 2*selene_c) RELATIVE arm_selene1
ROTATED (0.0,0.0,0.0) RELATIVE arm_selene1
COMPONENT PP_small_before_middle_focus = PSD_monitor(
filename = "PP_small_before_middle_focus",
nx = 200, ny = 200,
xwidth = 0.1, yheight = 0.1,
restore_neutron = 1)
AT (0, 0, 2*selene_c) RELATIVE arm_selene1
/* The proximal polariser comes next */
COMPONENT polarizer2 = Polariser(nIncRefr=1, d_substrate = 5e-4, reflect_d=0, reflect_u=0, lin=polarizer_start, length=polarizer_length,
delta_theta=(Theta1_polarizer+Theta2_polarizer)*PI/180.0, h2=0.1, h1=0.05, abs_ref=1, m_u=4.0, m_d=0.65, both_coated=1, alpha=2.3, W = 0.0014)
WHEN (enable_polarizer>1)
AT (0, 0, polarizer_start) RELATIVE arm_polarizer
ROTATED (0.0,0,90.0) RELATIVE arm_polarizer
/**************************
* Selene 2 neutron guide *
**************************/
/* Absorber to cut direct view beam (Copper in CAD model) */
COMPONENT slit_before_selene_guide_2 = Slit(
xmin = -selene_entry-0.005, xmax=0,
ymin = -selene_entry-0.01, ymax = 0)
AT (0, 0, selene_distance-0.002) RELATIVE arm_selene2
COMPONENT block_before_selene_guide_2 = Absorber(
xmin=-1, xmax=1,
ymin=-selene_entry/4.0+0.015,
ymax=selene_entry/4.0-0.015,
zmin=0.0, zmax=0.001)
AT (0, 0, selene_distance-0.0015) RELATIVE arm_selene2
/* Selene 2 elliptic guide first half */
COMPONENT selene_guide_21 = Elliptic_guide_gravity(
l=0.5*selene_length-0.001, dimensionsAt = "mid",
linyh = selene_distance, loutyh= selene_c+0.001,
linxw = selene_distance, loutxw= selene_c+0.001,
xwidth=selene_b*2, yheight=selene_b*2,
mright=selene_coating, mleft=0, mtop=0, mbottom=selene_coating,
enableGravity=enable_gravity)
AT (0, 0, selene_distance) RELATIVE arm_selene2
/* Absorber to cut direct view beam (Copper in CAD model) */
COMPONENT slit_within_selene_guide_2 = Slit(
xmin = -selene_b, xmax=-selene_b*0.45+0.005,
ymin = -selene_b, ymax=-selene_b*0.45+0.005)
AT (0, 0, 1.0*selene_c) RELATIVE arm_selene2
/*COMPONENT block_within_selene_guide_2 = Absorber(
xmin =-1, xmax=1,
ymin =-selene_b*0.25+0.001,
ymax = selene_b*0.25-0.001,
zmin=0.0, zmax=0.0001)
AT (0, 0, 3.0*selene_c+0.0001) RELATIVE arm_selene2*/
/* Selene 2 elliptic guide first half */
COMPONENT selene_guide_22 = Elliptic_guide_gravity(
l=0.5*selene_length-0.001, dimensionsAt = "mid",
linyh = selene_c+0.001, loutyh= selene_distance,
linxw = selene_c+0.001, loutxw= selene_distance,
xwidth=selene_b*2, yheight=selene_b*2,
mright=selene_coating, mleft=0, mtop=0, mbottom=selene_coating,
enableGravity=enable_gravity)
AT (0, 0, 1.0*selene_c+0.001) RELATIVE arm_selene2
/* Absorber to cut direct view beam (Copper in CAD model) */
COMPONENT slit_after_selene_guide_2 = Slit(
xmin = -selene_entry-0.001, xmax=0, xmax = 0.0,
ymin = -selene_entry-0.005, ymax = 0)
AT (0, 0, 2*selene_c-selene_distance+0.001) RELATIVE arm_selene2
COMPONENT block_after_selene_guide_2 = Absorber(
xmin=-1, xmax=1,
ymin=-selene_entry/4.0+0.015,
ymax=selene_entry/4.0-0.015,
zmin=0.0, zmax=0.001)
AT (0, 0, 2*selene_c-selene_distance+0.0015) RELATIVE arm_selene2
REMOVABLE COMPONENT FocusMonitor = PSD_monitor(
filename = "Focus",
nx = 100, ny = 100,
xwidth = 0.005, yheight = 0.005,
restore_neutron = 1)
AT (0, 0, 4*selene_c) RELATIVE ISCS
ROTATED (0, 1.25, 0) RELATIVE ISCS
FINALLY
%{
%}
END
+270
View File
@@ -0,0 +1,270 @@
/*******************************************************************************
* McStas instrument definition URL=http://www.mcstas.org
*
* Instrument: Estia_selene1
*
* %Identification
* Written by: Artur Glavic (artur.glavic@psi.ch); Jochen Stahn (jochen.stahn@psi.ch); Christine Klauser (christine.klauser@psi.ch)
* Date: 01. 03. 2018
* Origin: PSI
* Release: McStas 2.4.1
* Version: 1.0
* %INSTRUMENT_SITE: ESS (E02)
*
* Estia is a vertical sample, focusing reflectometer for small sample
*
* %Description
* This is the Selene guide 1 part of the instrument, it is not intended for use but is
* included in the Estia instrument model. The removable components, however,
* allow for a stand alone use for debugging purpose.
*
* %Parameters
* enable_gravity: [0/1] Use gravity in elliptical guide model.
*
* %End
*******************************************************************************/
DEFINE INSTRUMENT Estia_selene(int enable_gravity=0)
DECLARE
%{
// Selene 1
// Selene 2
/* selene parameters */
double selene_xi = 0.6 ; // length of neutron guide over 2*c of ellipses, xi
double selene_b = 0.1047; // m selene short half axis c
double selene_c = 6.0000; // m selene focus half distance c
/* Coating along the elliptic guide to allow sharp wavelenth cut-off */
double scoating_01 = 3.58; // coating at Selene 1 start/ Selene 2 end
double scoating_02 = 3.36;
double scoating_03 = 3.19;
double scoating_04 = 3.07;
double scoating_05 = 2.98;
double scoating_06 = 2.92;
double scoating_07 = 2.88;
double scoating_08 = 2.86; // coating in center
double scoating_09 = 2.86;
double scoating_10 = 2.87;
double scoating_11 = 2.91;
double scoating_12 = 2.98;
double scoating_13 = 3.07;
double scoating_14 = 3.20;
double scoating_15 = 3.38; // coating at Selene 1 end/ Selene 2 start
/* other variables */
double selene_entry ;
double selene_length ;
double selene_distance ;
double selene_segment ;
%}
INITIALIZE
%{
/* derived parameters for the selene guides */
selene_length = 2.0 * selene_c * selene_xi;
selene_distance = ( 1.0 - selene_xi ) * selene_c;
selene_entry = selene_b * sqrt( 1.0 - pow(selene_xi,2.0) );
selene_segment = selene_length / 15.;
%}
TRACE
REMOVABLE COMPONENT origin = Progress_bar()
AT (0,0,0) ABSOLUTE
REMOVABLE COMPONENT ISCS = Arm()
AT (0, 0, 0) RELATIVE origin
REMOVABLE COMPONENT arm_selene1 = Arm()
AT (0, 0, 0) RELATIVE ISCS
ROTATED (0, 0, 0) RELATIVE ISCS
REMOVABLE COMPONENT arm_selene2 = Arm()
AT (0, 0, 2*selene_c) RELATIVE ISCS
REMOVABLE COMPONENT source = Moderator(radius = 0.001, focus_xw = selene_entry+0.002, focus_yh = selene_entry+0.002,
target_index=3, Emin=1, Emax=15)
AT (0,0,0) RELATIVE ISCS
ROTATED (-1.25, 1.25, 0) RELATIVE ISCS
/**************************
* Selene 1 neutron guide *
**************************/
/* Absorber to cut direct view beam (Copper in CAD model) */
COMPONENT slit_before_selene_guide_1 = Slit(
xmin = 0, xmax = selene_entry+0.001,
ymin = -selene_entry-0.005, ymax = 0)
AT (0, 0, selene_distance-0.01) RELATIVE arm_selene1
COMPONENT block_before_selene_guide_1 = Absorber(
xmin=-1, xmax=1,
ymin=-selene_entry/4.0+0.005,
ymax=selene_entry/4.0-0.005,
zmin=0.0, zmax=0.001)
AT (0, 0, selene_distance-0.0095) RELATIVE arm_selene1
/* Selene 1 elliptic guide */
COMPONENT E02_03_01 = Elliptic_guide_gravity(
l=selene_segment, dimensionsAt = "mid",
linyh = selene_distance, loutyh= selene_distance+14*selene_segment,
linxw = selene_distance, loutxw= selene_distance+14*selene_segment,
xwidth=selene_b*2, yheight=selene_b*2,
mright=0, mleft=scoating_01, mbottom=scoating_01, mtop=0,
enableGravity=enable_gravity)
AT (0, 0, selene_distance) RELATIVE arm_selene1
COMPONENT E02_03_02 = Elliptic_guide_gravity(
l=selene_segment, dimensionsAt = "mid",
linyh = selene_distance+1*selene_segment, loutyh= selene_distance+13*selene_segment,
linxw = selene_distance+1*selene_segment, loutxw= selene_distance+13*selene_segment,
xwidth=selene_b*2, yheight=selene_b*2,
mright=0, mleft=scoating_02, mbottom=scoating_02, mtop=0,
enableGravity=enable_gravity)
AT (0, 0, selene_distance+1*selene_segment) RELATIVE arm_selene1
COMPONENT E02_03_03 = Elliptic_guide_gravity(
l=selene_segment, dimensionsAt = "mid",
linyh = selene_distance+2*selene_segment, loutyh= selene_distance+12*selene_segment,
linxw = selene_distance+2*selene_segment, loutxw= selene_distance+12*selene_segment,
xwidth=selene_b*2, yheight=selene_b*2,
mright=0, mleft=scoating_03, mbottom=scoating_03, mtop=0,
enableGravity=enable_gravity)
AT (0, 0, selene_distance+2*selene_segment) RELATIVE arm_selene1
COMPONENT E02_03_04 = Elliptic_guide_gravity(
l=selene_segment, dimensionsAt = "mid",
linyh = selene_distance+3*selene_segment, loutyh= selene_distance+11*selene_segment,
linxw = selene_distance+3*selene_segment, loutxw= selene_distance+11*selene_segment,
xwidth=selene_b*2, yheight=selene_b*2,
mright=0, mleft=scoating_04, mbottom=scoating_04, mtop=0,
enableGravity=enable_gravity)
AT (0, 0, selene_distance+3*selene_segment) RELATIVE arm_selene1
COMPONENT E02_03_05 = Elliptic_guide_gravity(
l=selene_segment, dimensionsAt = "mid",
linyh = selene_distance+4*selene_segment, loutyh= selene_distance+10*selene_segment,
linxw = selene_distance+4*selene_segment, loutxw= selene_distance+10*selene_segment,
xwidth=selene_b*2, yheight=selene_b*2,
mright=0, mleft=scoating_05, mbottom=scoating_05, mtop=0,
enableGravity=enable_gravity)
AT (0, 0, selene_distance+4*selene_segment) RELATIVE arm_selene1
COMPONENT E02_03_06 = Elliptic_guide_gravity(
l=selene_segment, dimensionsAt = "mid",
linyh = selene_distance+5*selene_segment, loutyh= selene_distance+9*selene_segment,
linxw = selene_distance+5*selene_segment, loutxw= selene_distance+9*selene_segment,
xwidth=selene_b*2, yheight=selene_b*2,
mright=0, mleft=scoating_06, mbottom=scoating_06, mtop=0,
enableGravity=enable_gravity)
AT (0, 0, selene_distance+5*selene_segment) RELATIVE arm_selene1
COMPONENT E02_03_07 = Elliptic_guide_gravity(
l=selene_segment, dimensionsAt = "mid",
linyh = selene_distance+6*selene_segment, loutyh= selene_distance+8*selene_segment,
linxw = selene_distance+6*selene_segment, loutxw= selene_distance+8*selene_segment,
xwidth=selene_b*2, yheight=selene_b*2,
mright=0, mleft=scoating_07, mbottom=scoating_07, mtop=0,
enableGravity=enable_gravity)
AT (0, 0, selene_distance+6*selene_segment) RELATIVE arm_selene1
COMPONENT E02_03_08 = Elliptic_guide_gravity(
l=selene_segment, dimensionsAt = "mid",
linyh = selene_distance+7*selene_segment, loutyh= selene_distance+7*selene_segment,
linxw = selene_distance+7*selene_segment, loutxw= selene_distance+7*selene_segment,
xwidth=selene_b*2, yheight=selene_b*2,
mright=0, mleft=scoating_08, mbottom=scoating_08, mtop=0,
enableGravity=enable_gravity)
AT (0, 0, selene_distance+7*selene_segment) RELATIVE arm_selene1
COMPONENT E02_03_09 = Elliptic_guide_gravity(
l=selene_segment, dimensionsAt = "mid",
linyh = selene_distance+8*selene_segment, loutyh= selene_distance+6*selene_segment,
linxw = selene_distance+8*selene_segment, loutxw= selene_distance+6*selene_segment,
xwidth=selene_b*2, yheight=selene_b*2,
mright=0, mleft=scoating_09, mbottom=scoating_09, mtop=0,
enableGravity=enable_gravity)
AT (0, 0, selene_distance+8*selene_segment) RELATIVE arm_selene1
COMPONENT E02_03_10 = Elliptic_guide_gravity(
l=selene_segment, dimensionsAt = "mid",
linyh = selene_distance+9*selene_segment, loutyh= selene_distance+5*selene_segment,
linxw = selene_distance+9*selene_segment, loutxw= selene_distance+5*selene_segment,
xwidth=selene_b*2, yheight=selene_b*2,
mright=0, mleft=scoating_10, mbottom=scoating_10, mtop=0,
enableGravity=enable_gravity)
AT (0, 0, selene_distance+9*selene_segment) RELATIVE arm_selene1
COMPONENT E02_03_11 = Elliptic_guide_gravity(
l=selene_segment, dimensionsAt = "mid",
linyh = selene_distance+10*selene_segment, loutyh= selene_distance+4*selene_segment,
linxw = selene_distance+10*selene_segment, loutxw= selene_distance+4*selene_segment,
xwidth=selene_b*2, yheight=selene_b*2,
mright=0, mleft=scoating_11, mbottom=scoating_11, mtop=0,
enableGravity=enable_gravity)
AT (0, 0, selene_distance+10*selene_segment) RELATIVE arm_selene1
COMPONENT E02_03_12 = Elliptic_guide_gravity(
l=selene_segment, dimensionsAt = "mid",
linyh = selene_distance+11*selene_segment, loutyh= selene_distance+3*selene_segment,
linxw = selene_distance+11*selene_segment, loutxw= selene_distance+3*selene_segment,
xwidth=selene_b*2, yheight=selene_b*2,
mright=0, mleft=scoating_12, mbottom=scoating_12, mtop=0,
enableGravity=enable_gravity)
AT (0, 0, selene_distance+11*selene_segment) RELATIVE arm_selene1
COMPONENT E02_03_13 = Elliptic_guide_gravity(
l=selene_segment, dimensionsAt = "mid",
linyh = selene_distance+12*selene_segment, loutyh= selene_distance+2*selene_segment,
linxw = selene_distance+12*selene_segment, loutxw= selene_distance+2*selene_segment,
xwidth=selene_b*2, yheight=selene_b*2,
mright=0, mleft=scoating_13, mbottom=scoating_13, mtop=0,
enableGravity=enable_gravity)
AT (0, 0, selene_distance+12*selene_segment) RELATIVE arm_selene1
COMPONENT E02_03_14 = Elliptic_guide_gravity(
l=selene_segment, dimensionsAt = "mid",
linyh = selene_distance+13*selene_segment, loutyh= selene_distance+1*selene_segment,
linxw = selene_distance+13*selene_segment, loutxw= selene_distance+1*selene_segment,
xwidth=selene_b*2, yheight=selene_b*2,
mright=0, mleft=scoating_14, mbottom=scoating_14, mtop=0,
enableGravity=enable_gravity)
AT (0, 0, selene_distance+13*selene_segment) RELATIVE arm_selene1
COMPONENT E02_03_15 = Elliptic_guide_gravity(
l=selene_segment, dimensionsAt = "mid",
linyh = selene_distance+14*selene_segment, loutyh= selene_distance,
linxw = selene_distance+14*selene_segment, loutxw= selene_distance,
xwidth=selene_b*2, yheight=selene_b*2,
mright=0, mleft=scoating_15, mbottom=scoating_15, mtop=0,
enableGravity=enable_gravity)
AT (0, 0, selene_distance+14*selene_segment) RELATIVE arm_selene1
/* Absorber to cut direct view beam (Copper in CAD model) */
COMPONENT slit_after_selene_guide_1 = Slit(
xmin = 0, xmax = selene_entry+0.005,
ymin = -selene_entry-0.005, ymax = 0)
AT (0, 0, 2*selene_c-selene_distance+0.001) RELATIVE arm_selene1
COMPONENT block_after_selene_guide_1 = Absorber(
xmin=-1, xmax=1,
ymin=-selene_entry/4.0+0.015,
ymax=selene_entry/4.0-0.015,
zmin=0.0, zmax=0.001)
AT (0, 0, 2*selene_c-selene_distance+0.0015) RELATIVE arm_selene1
FINALLY
%{
%}
END
+210
View File
@@ -0,0 +1,210 @@
/*******************************************************************************
* McStas instrument definition URL=http://www.mcstas.org
*
* Instrument: Estia_selene2
*
* %Identification
* Written by: Artur Glavic (artur.glavic@psi.ch); Jochen Stahn (jochen.stahn@psi.ch); Christine Klauser (christine.klauser@psi.ch)
* Date: 01. 03. 2018
* Origin: PSI
* Release: McStas 2.4.1
* Version: 1.0
* %INSTRUMENT_SITE: ESS (E02)
*
* Estia is a vertical sample, focusing reflectometer for small sample
*
* %Description
* This is the second Selene guide part of the instrument, it is not intended for use but is
* included in the Estia instrument model. This file does not work no its own, as it
* requires the definitions specified in the Estia_selene1.instr DECLAR and INITIALIZE.
*
* %End
*******************************************************************************/
DEFINE INSTRUMENT Estia_selene2()
DECLARE
%{
%}
INITIALIZE
%{
%}
TRACE
/**************************
* Selene 2 neutron guide *
**************************/
/* Absorber to cut direct view beam (Bor-Al in CAD model) */
COMPONENT slit_before_selene_guide_2 = Slit(
xmin = -selene_entry-0.005, xmax=0,
ymin = 0, ymax = selene_entry+0.01)
AT (0, 0, selene_distance-0.002) RELATIVE arm_selene2
COMPONENT block_before_selene_guide_2 = Absorber(
xmin=-1, xmax=1,
ymin=-selene_entry/4.0+0.015,
ymax=selene_entry/4.0-0.015,
zmin=0.0, zmax=0.001)
AT (0, 0, selene_distance-0.0015) RELATIVE arm_selene2
/* Selene 2 elliptic guide first half */
COMPONENT E02_04_01 = Elliptic_guide_gravity(
l=selene_segment, dimensionsAt = "mid",
linyh = selene_distance, loutyh= selene_distance+14*selene_segment,
linxw = selene_distance, loutxw= selene_distance+14*selene_segment,
xwidth=selene_b*2, yheight=selene_b*2,
mright=scoating_15, mleft=0, mbottom=0, mtop=scoating_15,
enableGravity=enable_gravity)
AT (0, 0, selene_distance) RELATIVE arm_selene2
COMPONENT E02_04_02 = Elliptic_guide_gravity(
l=selene_segment, dimensionsAt = "mid",
linyh = selene_distance+1*selene_segment, loutyh= selene_distance+13*selene_segment,
linxw = selene_distance+1*selene_segment, loutxw= selene_distance+13*selene_segment,
xwidth=selene_b*2, yheight=selene_b*2,
mright=scoating_14, mleft=0, mbottom=0, mtop=scoating_14,
enableGravity=enable_gravity)
AT (0, 0, selene_distance+1*selene_segment) RELATIVE arm_selene2
COMPONENT E02_04_03 = Elliptic_guide_gravity(
l=selene_segment, dimensionsAt = "mid",
linyh = selene_distance+2*selene_segment, loutyh= selene_distance+12*selene_segment,
linxw = selene_distance+2*selene_segment, loutxw= selene_distance+12*selene_segment,
xwidth=selene_b*2, yheight=selene_b*2,
mright=scoating_13, mleft=0, mbottom=0, mtop=scoating_13,
enableGravity=enable_gravity)
AT (0, 0, selene_distance+2*selene_segment) RELATIVE arm_selene2
COMPONENT E02_04_04 = Elliptic_guide_gravity(
l=selene_segment, dimensionsAt = "mid",
linyh = selene_distance+3*selene_segment, loutyh= selene_distance+11*selene_segment,
linxw = selene_distance+3*selene_segment, loutxw= selene_distance+11*selene_segment,
xwidth=selene_b*2, yheight=selene_b*2,
mright=scoating_12, mleft=0, mbottom=0, mtop=scoating_12,
enableGravity=enable_gravity)
AT (0, 0, selene_distance+3*selene_segment) RELATIVE arm_selene2
COMPONENT E02_04_05 = Elliptic_guide_gravity(
l=selene_segment, dimensionsAt = "mid",
linyh = selene_distance+4*selene_segment, loutyh= selene_distance+10*selene_segment,
linxw = selene_distance+4*selene_segment, loutxw= selene_distance+10*selene_segment,
xwidth=selene_b*2, yheight=selene_b*2,
mright=scoating_11, mleft=0, mbottom=0, mtop=scoating_11,
enableGravity=enable_gravity)
AT (0, 0, selene_distance+4*selene_segment) RELATIVE arm_selene2
COMPONENT E02_04_06 = Elliptic_guide_gravity(
l=selene_segment, dimensionsAt = "mid",
linyh = selene_distance+5*selene_segment, loutyh= selene_distance+9*selene_segment,
linxw = selene_distance+5*selene_segment, loutxw= selene_distance+9*selene_segment,
xwidth=selene_b*2, yheight=selene_b*2,
mright=scoating_10, mleft=0, mbottom=0, mtop=scoating_10,
enableGravity=enable_gravity)
AT (0, 0, selene_distance+5*selene_segment) RELATIVE arm_selene2
COMPONENT E02_04_07 = Elliptic_guide_gravity(
l=selene_segment, dimensionsAt = "mid",
linyh = selene_distance+6*selene_segment, loutyh= selene_distance+8*selene_segment,
linxw = selene_distance+6*selene_segment, loutxw= selene_distance+8*selene_segment,
xwidth=selene_b*2, yheight=selene_b*2,
mright=scoating_09, mleft=0, mbottom=0, mtop=scoating_09,
enableGravity=enable_gravity)
AT (0, 0, selene_distance+6*selene_segment) RELATIVE arm_selene2
COMPONENT E02_04_08 = Elliptic_guide_gravity(
l=selene_segment, dimensionsAt = "mid",
linyh = selene_distance+7*selene_segment, loutyh= selene_distance+7*selene_segment,
linxw = selene_distance+7*selene_segment, loutxw= selene_distance+7*selene_segment,
xwidth=selene_b*2, yheight=selene_b*2,
mright=scoating_08, mleft=0, mbottom=0, mtop=scoating_08,
enableGravity=enable_gravity)
AT (0, 0, selene_distance+7*selene_segment) RELATIVE arm_selene2
COMPONENT E02_04_09 = Elliptic_guide_gravity(
l=selene_segment, dimensionsAt = "mid",
linyh = selene_distance+8*selene_segment, loutyh= selene_distance+6*selene_segment,
linxw = selene_distance+8*selene_segment, loutxw= selene_distance+6*selene_segment,
xwidth=selene_b*2, yheight=selene_b*2,
mright=scoating_07, mleft=0, mbottom=0, mtop=scoating_07,
enableGravity=enable_gravity)
AT (0, 0, selene_distance+8*selene_segment) RELATIVE arm_selene2
COMPONENT E02_04_10 = Elliptic_guide_gravity(
l=selene_segment, dimensionsAt = "mid",
linyh = selene_distance+9*selene_segment, loutyh= selene_distance+5*selene_segment,
linxw = selene_distance+9*selene_segment, loutxw= selene_distance+5*selene_segment,
xwidth=selene_b*2, yheight=selene_b*2,
mright=scoating_06, mleft=0, mbottom=0, mtop=scoating_06,
enableGravity=enable_gravity)
AT (0, 0, selene_distance+9*selene_segment) RELATIVE arm_selene2
COMPONENT E02_04_11 = Elliptic_guide_gravity(
l=selene_segment, dimensionsAt = "mid",
linyh = selene_distance+10*selene_segment, loutyh= selene_distance+4*selene_segment,
linxw = selene_distance+10*selene_segment, loutxw= selene_distance+4*selene_segment,
xwidth=selene_b*2, yheight=selene_b*2,
mright=scoating_05, mleft=0, mbottom=0, mtop=scoating_05,
enableGravity=enable_gravity)
AT (0, 0, selene_distance+10*selene_segment) RELATIVE arm_selene2
COMPONENT E02_04_12 = Elliptic_guide_gravity(
l=selene_segment, dimensionsAt = "mid",
linyh = selene_distance+11*selene_segment, loutyh= selene_distance+3*selene_segment,
linxw = selene_distance+11*selene_segment, loutxw= selene_distance+3*selene_segment,
xwidth=selene_b*2, yheight=selene_b*2,
mright=scoating_04, mleft=0, mbottom=0, mtop=scoating_04,
enableGravity=enable_gravity)
AT (0, 0, selene_distance+11*selene_segment) RELATIVE arm_selene2
COMPONENT E02_04_13 = Elliptic_guide_gravity(
l=selene_segment, dimensionsAt = "mid",
linyh = selene_distance+12*selene_segment, loutyh= selene_distance+2*selene_segment,
linxw = selene_distance+12*selene_segment, loutxw= selene_distance+2*selene_segment,
xwidth=selene_b*2, yheight=selene_b*2,
mright=scoating_03, mleft=0, mbottom=0, mtop=scoating_03,
enableGravity=enable_gravity)
AT (0, 0, selene_distance+12*selene_segment) RELATIVE arm_selene2
COMPONENT E02_04_14 = Elliptic_guide_gravity(
l=selene_segment, dimensionsAt = "mid",
linyh = selene_distance+13*selene_segment, loutyh= selene_distance+1*selene_segment,
linxw = selene_distance+13*selene_segment, loutxw= selene_distance+1*selene_segment,
xwidth=selene_b*2, yheight=selene_b*2,
mright=scoating_02, mleft=0, mbottom=0, mtop=scoating_02,
enableGravity=enable_gravity)
AT (0, 0, selene_distance+13*selene_segment) RELATIVE arm_selene2
COMPONENT E02_04_15 = Elliptic_guide_gravity(
l=selene_segment, dimensionsAt = "mid",
linyh = selene_distance+14*selene_segment, loutyh= selene_distance,
linxw = selene_distance+14*selene_segment, loutxw= selene_distance,
xwidth=selene_b*2, yheight=selene_b*2,
mright=scoating_01, mleft=0, mbottom=0, mtop=scoating_01,
enableGravity=enable_gravity)
AT (0, 0, selene_distance+14*selene_segment) RELATIVE arm_selene2
/* Absorber to cut direct view beam (Bor-Al in CAD model) */
COMPONENT slit_after_selene_guide_2 = Slit(
xmin = -selene_entry-0.001, xmax=0, xmax = 0.0,
ymin = 0, ymax = selene_entry+0.01)
AT (0, 0, 2*selene_c-selene_distance+0.001) RELATIVE arm_selene2
COMPONENT block_after_selene_guide_2 = Absorber(
xmin=-1, xmax=1,
ymin=-selene_entry/4.0+0.015,
ymax=selene_entry/4.0-0.015,
zmin=0.0, zmax=0.001)
AT (0, 0, 2*selene_c-selene_distance+0.0015) RELATIVE arm_selene2
FINALLY
%{
%}
END
+91
View File
@@ -0,0 +1,91 @@
/*******************************************************************************
*
* McStas, neutron ray-tracing package
* Copyright 1997-2002, All rights reserved
* Risoe National Laboratory, Roskilde, Denmark
* Institut Laue Langevin, Grenoble, France
*
* Written by: Erik B Knudsen
* Date: May 2017
* Version: $Revision$
* Release: McStas 2.4
* Origin: DTU Physics
*
* Component: Pol_RFSF_ideal
*
*
* %I
*
* Ideal model of a spin flipper
*
* %D
* This component simply mirrors the polarization vector of the neutron
* ray in the plane through (0,0,0) with normal nx,ny,nz.
* The flipper is surrounded by a perfectly absorbing box. Neutron rays not hitting
* the box are left untouched.
*
* %P
* Input parameters:
* nx: [ ] x-component of the normal vector of the flipping plane.
* ny: [ ] y-component of the normal vector of the flipping plane.
* nz: [ ] z-component of the normal vector of the flipping plane.
* xwidth: [m] width of the spin flipper.
* yheight: [m] height of the spin flipper.
* zdepth: [m] length of the spin flipper.
*
*
* %E
*******************************************************************************/
DEFINE COMPONENT Pol_SF_ideal
DEFINITION PARAMETERS ()
SETTING PARAMETERS (nx=0, ny=1, nz=0, xwidth=0.1, yheight=0.1, zdepth=0.1)
OUTPUT PARAMETERS ()
/* Neutron parameters: (x,y,z,vx,vy,vz,t,sx,sy,sz,p) */
TRACE
%{
int hit;
double t0,t1;
hit=box_intersect(&t0,&t1, x,y,z,vx,vy,vz, xwidth,yheight,zdepth);
if(hit){
PROP_DT(t0);
if(fabs(z- -zdepth*0.5)>DBL_EPSILON){
/*neutron must have hit the side walls*/
ABSORB;
}
/*move to center of box and flip*/
PROP_Z0;
SCATTER;
double s=scalar_prod(sx,sy,sz,nx,ny,nz);
if (s!=0){
sx-=2*s*nx;
sy-=2*s*ny;
sz-=2*s*nz;
}
PROP_DT((t1-t0)/2);/*propagate the remaining distance to the box exit*/
if(fabs(z-zdepth*0.5)>DBL_EPSILON){
ABSORB;
}
}
%}
MCDISPLAY
%{
double dx=xwidth/16;
double dy=yheight/8;
box(0,0,0,xwidth,yheight,zdepth);
line(dx,-dy,0,dx,-dy+yheight/2.0,0);
line(-dx,dy,0,-dx,dy-yheight/2.0,0);
line(dx,-dy+yheight/2.0,0, dx+xwidth/16,-dy+yheight-yheight/16,0);
line(dx,-dy+yheight/2.0,0, dx-xwidth/16,-dy+yheight-yheight/16,0);
line(-dx,dy-yheight/2.0,0, -dx+xwidth/16,dy-yheight+yheight/16,0);
line(-dx,dy-yheight/2.0,0, -dx-xwidth/16,dy-yheight+yheight/16,0);
%}
END
+605 -488
View File
File diff suppressed because it is too large Load Diff
+108
View File
@@ -0,0 +1,108 @@
/*******************************************************************************
* McStas instrument definition URL=http://www.mcstas.org
*
* Instrument: PolarizerTest
*
* %Identification
* Written by: Artur Glavic (artur.glavic@psi.ch)
* Date: 25. 07. 2019
* Origin: PSI
* Release: McStas 2.5
* Version: 1.0
*
* Test for polarization parameters of curved transmission polarizer
*
* %Description
*
*
* %Parameters
*
* %End
*******************************************************************************/
DEFINE INSTRUMENT PolarizerTest(m_up=5.08, m_down=0.6, R0=1.00,
m_residual=0.55, Theta=1.5,
lambda_min=2.0, lambda_max=32.0,
gamma=1.66, focus_height=0.00001,
phi=0.0, polarizer_start=0.2
)
DECLARE
%{
double lambda_mean,delta_lambda,polarizer_length,zero_position;
double source_div, detector_height;
int use_streight=0;
%}
INITIALIZE
%{
lambda_mean=(lambda_min+lambda_max)/2.0;
delta_lambda=(lambda_max-lambda_min)/2.0;
polarizer_length=(exp(Theta/gamma)-1.0)*polarizer_start;
zero_position=exp(0.5*Theta/gamma)*polarizer_start;
source_div=0.1*Theta/180.*PI;
detector_height=source_div/0.1*(polarizer_length+polarizer_start+0.5);
%}
TRACE
COMPONENT origin = Progress_bar()
AT (0,0,0) ABSOLUTE
COMPONENT focus = Arm()
AT (0,0,0) RELATIVE origin
ROTATED (phi,0,0) RELATIVE origin
COMPONENT Source = Source_simple(xwidth = 0.00001, yheight = focus_height,
dist = 0.1, focus_xw = 0.0001,
focus_yh = source_div,
lambda0 = lambda_mean, dlambda = delta_lambda,
gauss=0)
AT (0, 0, 0) RELATIVE origin
COMPONENT FluxIn = L_monitor(xwidth = 0.5, yheight = 0.5, filename="fluxin",
nL=601, Lmin=lambda_min, Lmax=lambda_max)
AT (0, 0, 0.101) RELATIVE origin
COMPONENT polarizer = Polariser(enable_ref=1, abs_ref=0,
reflect_d=0, reflect_u=0, T_loss=4.0e3,
m_substrate=0.469522, d_substrate = 5e-4,
m_u=m_up, m_d=m_down, R0=R0,
m_residual=m_residual,
lin=polarizer_start, length=polarizer_length,
delta_theta=(Theta)*PI/180.0, h2=0.1, h1=0.05,
both_coated=1, alpha=2.3, W = 0.0014)
WHEN use_streight==0
AT (0, -5e-4, polarizer_start) RELATIVE focus
COMPONENT streight = Mirror(xwidth = 0.2, yheight = 0.2,
m=5.0, center=1, transmit=0)
WHEN use_streight==1
AT (0, 0, zero_position) RELATIVE origin
ROTATED (90.0-gamma,0,0) RELATIVE origin
COMPONENT FluxOut = L_monitor(xwidth = 0.5, yheight = detector_height,
filename="fluxout",
nL=601, Lmin=lambda_min, Lmax=lambda_max)
AT (0, 0, polarizer_length+polarizer_start+0.5) RELATIVE origin
COMPONENT PolMon = MeanPolLambda_monitor(xwidth=0.5, yheight=detector_height,
mx=-1,
nL=601, Lmin=lambda_min, Lmax=lambda_max,
filename="polarization")
AT (0, 0, polarizer_length+polarizer_start+0.5) RELATIVE origin
FINALLY
%{
%}
END
+58
View File
@@ -0,0 +1,58 @@
4.37173E-05 0.253838498
0.001967279 0.548126562
0.004240578 0.781912929
0.006513878 0.831781609
0.008656026 0.935568791
0.010841891 0.996046699
0.012984039 0.994254859
0.015038752 0.780745054
0.017443204 0.571166135
0.019497917 0.421992388
0.021771217 0.347878292
0.023913365 0.314863976
0.026317816 0.273648917
0.028241378 0.248910367
0.030514678 0.229613667
0.032700543 0.210315136
0.034711539 0.201079179
0.03711599 0.179153501
0.039301855 0.174572119
0.041444003 0.161570874
0.043629868 0.159330731
0.045772016 0.158627627
0.047957881 0.133760075
0.050231181 0.128230852
0.052373329 0.135460425
0.054559194 0.125901964
0.056701342 0.120064418
0.058974642 0.122922579
0.061160507 0.127084935
0.063171503 0.12641248
0.06548852 0.117637177
0.067630668 0.105706827
0.069903968 0.112415362
0.072089833 0.119429932
0.074231981 0.105676178
0.076417846 0.107716872
0.078428842 0.107129046
0.080833293 0.10262909
0.082888007 0.094164854
0.085161306 0.092239896
0.087303454 0.085655788
0.089620471 0.081615114
0.091631467 0.076984785
0.093817332 0.076157254
0.096090632 0.070972035
0.09823278 0.060551277
0.100287493 0.071366358
0.102691945 0.065426801
0.104834093 0.059045344
0.107019958 0.056229269
0.109162106 0.061290519
0.111347971 0.051160701
0.113490119 0.053836338
0.115763418 0.05839439
0.117949284 0.057514252
0.120091431 0.047536765
0.122364731 0.04877829
0.300000000 0.000000000
+58
View File
@@ -0,0 +1,58 @@
4.37173E-05 0.224586133
0.001967279 0.518057252
0.004240578 0.775672976
0.006513878 0.825879049
0.008656026 0.899985882
0.010841891 0.967789786
0.012984039 0.99196296
0.015038752 0.992909754
0.017443204 0.994177124
0.019497917 0.991424846
0.021771217 0.993887344
0.023913365 0.998648926
0.026317816 0.999361915
0.028241378 0.998678679
0.030514678 0.998109885
0.032700543 0.997772638
0.034711539 0.997392554
0.03711599 0.999600524
0.039301855 0.994779551
0.041444003 0.996274647
0.043629868 0.996149207
0.045772016 0.998305774
0.047957881 0.99658353
0.050231181 0.995950077
0.052373329 0.995563046
0.054559194 0.995479284
0.056701342 0.993347667
0.058974642 0.990508095
0.061160507 0.990476431
0.063171503 0.988754593
0.06548852 0.987959462
0.067630668 0.987058572
0.069903968 0.988939946
0.072089833 0.986821236
0.074231981 0.985149296
0.076417846 0.987362592
0.078428842 0.986507086
0.080833293 0.98295765
0.082888007 0.983313364
0.085161306 0.983817788
0.087303454 0.983554288
0.089620471 0.979763465
0.091631467 0.975348682
0.093817332 0.975887375
0.096090632 0.972236927
0.09823278 0.965515882
0.100287493 0.962852402
0.102691945 0.956266464
0.104834093 0.953509902
0.107019958 0.945236222
0.109162106 0.808306574
0.111347971 0.186001557
0.113490119 0.078929784
0.115763418 0.051976893
0.117949284 0.05338487
0.120091431 0.04844847
0.122364731 0.053399168
0.300000000 0.000000000
Binary file not shown.
+7 -3
View File
@@ -1,9 +1,13 @@
#!/bin/bash
export LD_LIBRARY_PATH=$LD_LIBRARY_PATH:/afs/psi.ch/project/sinq/sl6-64/mcstas2.4/mcstas/2.4/libs
if [ Estia_baseline.instr -nt Estia_baseline.out ] || [ ! -f Estia_baseline.out ] \
|| [ Estia_feeder.instr -nt Estia_baseline.out ] \
|| [ Estia_selene.instr -nt Estia_baseline.out ]; then
|| [ Estia_selene1.instr -nt Estia_baseline.out ] \
|| [ Estia_mf.instr -nt Estia_baseline.out ] \
|| [ Estia_selene2.instr -nt Estia_baseline.out ]; then
rm Estia_baseline.c Estia_baseline.out
mcstas -o Estia_baseline.c Estia_baseline.instr
mpicc -O3 -o Estia_baseline.out Estia_baseline.c -lm -DUSE_MPI -DUSE_NEXUS -lNeXus
mcstas -o Estia_baseline.c Estia_baseline.instr --trace
mpicc -O2 -o Estia_baseline.out Estia_baseline.c -lm -DUSE_MPI -DUSE_NEXUS -lNeXus
fi
+3
View File
@@ -0,0 +1,3 @@
mcc05:72
mcc06:72
+83
View File
@@ -0,0 +1,83 @@
#!/usr/bin/env python
'''
Run a simulation through various values of polarizer 1/2 distance to focus
and generate a datafile with analyzed values for polarization and transmission
for the short and long wavelength area.
'''
import sys, os
from numpy import *
sys.path.append('../analysis')
from mcstas_reader import McSim
Iref=1.0
L=linspace(1.95, 32.05, 302)
Lc=(L[:-1]+L[1:])/2.
gamma1=1.8
gamma2=1.66
command_s_ref="mpiexec -np 60 --hostfile hostnames Estia_baseline.out --ncount=1e10 --dir=../results/polarizer_s_ref --format=NeXus sample=1 omegaa=2.0 over_illumination=0.001 lambda_start=2.0 lambda_end=32.0 sample_length=0.01 enable_polarizer=0"
command_l_ref="mpiexec -np 60 --hostfile hostnames Estia_baseline.out --ncount=2e9 --dir=../results/polarizer_l_ref --format=NeXus sample=1 omegaa=6.0 over_illumination=0.001 lambda_start=2.0 lambda_end=32.0 sample_length=0.048 enable_polarizer=0"
command_s="mpiexec -np 60 --hostfile hostnames Estia_baseline.out --ncount=5e9 --dir=../results/polarizer_s --format=NeXus sample=1 omegaa=2.0 over_illumination=0.001 pol1_start=%%f pol1_angle=%f pol2_start=%%f pol2_angle=%f lambda_start=2.0 lambda_end=32.0 sample_length=0.01 enable_polarizer=3"%(gamma1, gamma2)
command_l="mpiexec -np 60 --hostfile hostnames Estia_baseline.out --ncount=1e9 --dir=../results/polarizer_l --format=NeXus sample=1 omegaa=6.0 over_illumination=0.001 pol1_start=%%f pol1_angle=%f pol2_start=%%f pol2_angle=%f lambda_start=2.0 lambda_end=32.0 sample_length=0.048 enable_polarizer=3"%(gamma1, gamma2)
lengths=linspace(0.1, 0.6, 25)
def get_polarization(ds):
N,ignore=histogram(ds.data['L'], bins=L)
PN,ignore=histogram(ds.data['L'], bins=L, weights=sqrt(ds.data['sx']**2+ds.data['sy']**2+ds.data['sz']**2))
return PN/maximum(N,1)
def get_intensity(ds):
I, ignore=histogram(ds.data['L'], bins=L, weights=ds.data['p'])
return I
def analyze(ds):
P=get_polarization(ds)
I=get_intensity(ds)
T=I/maximum(Iref, Iref[Iref>0].min())
return P,T
def get_data(fname):
det=McSim(fname)['tof_detector']
P,T=analyze(det)
return [P[17], P[230], T[17], T[230]]
def simulate_and_analyze(pol1_start, pol2_start):
print pol1_start, pol2_start
os.system('rm -rf ../results/polarizer_?')
os.system(command_s%(pol1_start, pol2_start))
os.system(command_l%(pol1_start, pol2_start))
global Iref
Iref=Iref_s
ds=get_data('../results/polarizer_s')
Iref=Iref_l
dl=get_data('../results/polarizer_l')
with open('../results/polarizer_data.dat', 'a') as fh:
fh.write((10*'%.6f '+'\n')%tuple([pol1_start, pol2_start]+ds+dl))
if __name__=='__main__':
if not (os.path.exists('../results/polarizer_s_ref') and
os.path.exists('../results/polarizer_l_ref')):
os.system('rm -rf ../results/polarizer_?_ref')
print 'Simulate reference intensities'
os.system(command_s_ref)
os.system(command_l_ref)
ref=McSim('../results/polarizer_s_ref')['tof_detector']
Iref_s=get_intensity(ref)
ref=McSim('../results/polarizer_l_ref')['tof_detector']
Iref_l=get_intensity(ref)
with open('../results/polarizer_data.dat', 'w') as fh:
fh.write('# pol1 pol2 P3.7_s P25_s T3.7_s T25_s P3.7_l P25_l T3.7_l T25_l\n')
for pol1 in lengths:
for pol2 in lengths:
simulate_and_analyze(pol1, pol2)
with open('../results/polarizer_data.dat', 'a') as fh:
fh.write('\n')
+7
View File
@@ -0,0 +1,7 @@
!#/bin/sh
rm PolarizerTest.c PolarizerTest.out
rm -r ../results/polarizer_coating
mcstas -o PolarizerTest.c PolarizerTest.instr --trace
#mpicc -O2 -o PolarizerTest.out PolarizerTest.c -lm -DUSE_MPI
mpicc -O2 -o PolarizerTest.out PolarizerTest.c -lm -DUSE_MPI -DUSE_NEXUS -lNeXus -I/usr/local/include/nexus
LD_LIBRARY_PATH=/usr/local/lib/ mpirun -np 6 ./PolarizerTest.out --ncount=1e8 --format=NeXus --dir=../results/polarizer_coating R0=1.0 gamma=1.66
+12
View File
@@ -0,0 +1,12 @@
#!/bin/bash
mcstas -o Estia_monitor.c Estia_monitor.instr
mpicc -O3 -o Estia_monitor.out Estia_monitor.c -lm -DUSE_MPI
mpirun -np 6 Estia_monitor.out --ncount=1e8 --dir=../results/monitor --gravitation enable_windows=1 lambda_start=1.0 lambda_end=35.0 direct_beam=0
mpirun -np 6 Estia_monitor.out --ncount=1e8 --dir=../results/monitor_ref --gravitation enable_windows=1 lambda_start=1.0 lambda_end=35.0 direct_beam=2
mpirun -np 6 Estia_monitor.out --ncount=1e8 --dir=../results/monitor_trans --gravitation enable_windows=1 lambda_start=1.0 lambda_end=35.0 direct_beam=1
mpirun -np 6 Estia_monitor.out --ncount=1e8 --dir=../results/monitor100 --gravitation enable_windows=1 lambda_start=1.0 lambda_end=35.0 direct_beam=0 foil_thickness=0.0001
mpirun -np 6 Estia_monitor.out --ncount=1e8 --dir=../results/monitor001 --gravitation enable_windows=1 lambda_start=1.0 lambda_end=35.0 direct_beam=0 foil_thickness=0.000001
+54 -190
View File
@@ -1,213 +1,77 @@
#!/bin/bash
DEST=../results
ncount=2e7
use_cores=4
ncount=1e8
use_cores=16
sample=4
omega=0.8
sample_length=0.01
sample_length=0.005
sample_height=0.01
lambda_start=3.0
lambda_end=32.0
###################### Brilliance Transfer 10x10mm² and 1x1mm² VS ####################
if [ Estia_baseline.instr -nt Estia_baseline.out ]; then
rm Estia_baseline.c Estia_baseline.out
mcstas -o Estia_baseline.c Estia_baseline.instr
mpicc -O3 -o Estia_baseline.out Estia_baseline.c -lm -DUSE_MPI -DUSE_NEXUS -lNeXus
fi
#
# if [ Estia_baseline_ana1.instr -nt Estia_baseline_ana1.out ]; then
# rm Estia_baseline_ana1.c Estia_baseline_ana1.out
# mcstas -o Estia_baseline_ana1.c Estia_baseline_ana1.instr
# mpicc -O3 -o Estia_baseline_ana1.out Estia_baseline_ana1.c -lm -DUSE_MPI -DUSE_NEXUS -lNeXus
# fi
#
# if [ Estia_baseline_ana2.instr -nt Estia_baseline_ana2.out ]; then
# rm Estia_baseline_ana2.c Estia_baseline_ana2.out
# mcstas -o Estia_baseline_ana2.c Estia_baseline_ana2.instr
# mpicc -O3 -o Estia_baseline_ana2.out Estia_baseline_ana2.c -lm -DUSE_MPI -DUSE_NEXUS -lNeXus
# fi
bash compile_if_needed.sh
###################### Reference measurement 10x5mm² and 10x0.5mm² VS ####################
###################### Reference and Ni-layer measurement 10x10mm² sample ####################
# ncount=1e10
sample_length=0.01
ncount=1e8
sample_length=0.005
sample_height=0.01
lambda_start=3.5
lambda_end=30.0
sample=1
# omega=1.0
# DESTi=$DEST/pol_ref_10x10_10
# if [ -e "$DESTi" ]; then
# rm -r "$DESTi"
# fi
#
# mpirun -np $use_cores Estia_baseline.out \
# --dir="$DESTi" --ncount=$ncount --gravitation \
# omegaa=$omega operationmode=0 theta_resolution=0.04 over_illumination=0.000 \
# sample=$sample sample_length=$sample_length sample_height=$sample_height \
# lambda_start=$lambda_start lambda_end=$lambda_end enable_gravity=1 enable_chopper=0 enable_polarizer=0 enable_analyzer=0
#
# DESTi=$DEST/pol1_10x10_10
# if [ -e "$DESTi" ]; then
# rm -r "$DESTi"
# fi
#
# mpirun -np $use_cores Estia_baseline.out \
# --dir="$DESTi" --ncount=$ncount --gravitation \
# omegaa=$omega operationmode=0 theta_resolution=0.04 over_illumination=0.000 \
# sample=$sample sample_length=$sample_length sample_height=$sample_height \
# lambda_start=$lambda_start lambda_end=$lambda_end enable_gravity=1 enable_chopper=0 enable_polarizer=1 enable_analyzer=0
#
# DESTi=$DEST/pol2_10x10_10
# if [ -e "$DESTi" ]; then
# rm -r "$DESTi"
# fi
#
# mpirun -np $use_cores Estia_baseline.out \
# --dir="$DESTi" --ncount=$ncount --gravitation \
# omegaa=$omega operationmode=0 theta_resolution=0.04 over_illumination=0.000 \
# sample=$sample sample_length=$sample_length sample_height=$sample_height \
# lambda_start=$lambda_start lambda_end=$lambda_end enable_gravity=1 enable_chopper=0 enable_polarizer=2 enable_analyzer=0
# #
# omega=7.0
# DESTi=$DEST/pol1_10x10_70
# if [ -e "$DESTi" ]; then
# rm -r "$DESTi"
# fi
#
# mpirun -np $use_cores Estia_baseline.out \
# --dir="$DESTi" --ncount=$ncount --gravitation \
# omegaa=$omega operationmode=0 theta_resolution=0.04 over_illumination=0.000 \
# sample=$sample sample_length=$sample_length sample_height=$sample_height \
# lambda_start=$lambda_start lambda_end=$lambda_end enable_gravity=1 enable_chopper=0 enable_polarizer=1 enable_analyzer=0
#
# DESTi=$DEST/pol2_10x10_70
# if [ -e "$DESTi" ]; then
# rm -r "$DESTi"
# fi
#
# mpirun -np $use_cores Estia_baseline.out \
# --dir="$DESTi" --ncount=$ncount --gravitation \
# omegaa=$omega operationmode=0 theta_resolution=0.04 over_illumination=0.000 \
# sample=$sample sample_length=$sample_length sample_height=$sample_height \
# lambda_start=$lambda_start lambda_end=$lambda_end enable_gravity=1 enable_chopper=0 enable_polarizer=2 enable_analyzer=0
#
# ncount=1e8
omega=0.0
# DESTi=$DEST/pol1_10x50_70
# if [ -e "$DESTi" ]; then
# rm -r "$DESTi"
# fi
#
# mpirun -np $use_cores Estia_baseline.out \
# --dir="$DESTi" --ncount=$ncount --gravitation \
# omegaa=$omega operationmode=0 theta_resolution=0.04 over_illumination=0.0025 \
# sample=$sample sample_length=$sample_length sample_height=$sample_height \
# lambda_start=$lambda_start lambda_end=$lambda_end enable_gravity=1 enable_chopper=0 enable_polarizer=1 enable_analyzer=0
DESTi=$DEST/pol2_10x50_70
DESTi=$DEST/pol_ref_10x5
if [ -e "$DESTi" ]; then
rm -r "$DESTi"
fi
mpirun -np $use_cores Estia_baseline.out \
--dir="$DESTi" --ncount=$ncount --gravitation \
omegaa=$omega operationmode=0 theta_resolution=0.04 over_illumination=0.0025 \
mpirun -np $use_cores --hostfile hostnames Estia_baseline.out \
--dir="$DESTi" --format=NeXus --ncount=$ncount \
sample=$sample sample_length=$sample_length sample_height=$sample_height \
lambda_start=$lambda_start lambda_end=$lambda_end enable_gravity=1 enable_chopper=0 enable_polarizer=2 enable_analyzer=0
#
# ncount=5e9
# sample_length=0.01
# omega=1.0
# sample_length=0.003
# DESTi=$DEST/pol_ref_3x10_10
# if [ -e "$DESTi" ]; then
# rm -r "$DESTi"
# fi
#
# mpirun -np $use_cores Estia_baseline.out \
# --dir="$DESTi" --ncount=$ncount --gravitation \
# omegaa=$omega operationmode=0 theta_resolution=0.04 over_illumination=0.000 \
# sample=$sample sample_length=$sample_length sample_height=$sample_height \
# lambda_start=$lambda_start lambda_end=$lambda_end enable_gravity=1 enable_chopper=0 enable_polarizer=0 enable_analyzer=0
#
# DESTi=$DEST/pol1_3x10_10
# if [ -e "$DESTi" ]; then
# rm -r "$DESTi"
# fi
#
# mpirun -np $use_cores Estia_baseline.out \
# --dir="$DESTi" --ncount=$ncount --gravitation \
# omegaa=$omega operationmode=0 theta_resolution=0.04 over_illumination=0.000 \
# sample=$sample sample_length=$sample_length sample_height=$sample_height \
# lambda_start=$lambda_start lambda_end=$lambda_end enable_gravity=1 enable_chopper=0 enable_polarizer=1 enable_analyzer=0
#
# DESTi=$DEST/pol2_3x10_10
# if [ -e "$DESTi" ]; then
# rm -r "$DESTi"
# fi
#
# mpirun -np $use_cores Estia_baseline.out \
# --dir="$DESTi" --ncount=$ncount --gravitation \
# omegaa=$omega operationmode=0 theta_resolution=0.04 over_illumination=0.000 \
# sample=$sample sample_length=$sample_length sample_height=$sample_height \
# lambda_start=$lambda_start lambda_end=$lambda_end enable_gravity=1 enable_chopper=0 enable_polarizer=2 enable_analyzer=0
#
#
# sample_height=0.01
# sample_length=0.01
# DESTi=$DEST/ana1_10x10_10
# if [ -e "$DESTi" ]; then
# rm -r "$DESTi"
# fi
#
# mpirun -np $use_cores Estia_baseline.out \
# --dir="$DESTi" --ncount=$ncount --gravitation \
# omegaa=$omega operationmode=0 theta_resolution=0.04 over_illumination=0.000 \
# sample=$sample sample_length=$sample_length sample_height=$sample_height \
# lambda_start=$lambda_start lambda_end=$lambda_end enable_gravity=1 enable_chopper=0 enable_polarizer=0 enable_analyzer=1
lambda_start=$lambda_start lambda_end=$lambda_end enable_gravity=0 enable_chopper=0 \
enable_polarizer=0 enable_analyzer=0
ncount=5e8
sample_length=0.0005
sample_height=0.01
# sample_height=0.001
# DESTi=$DEST/ana1_10x1_10
# if [ -e "$DESTi" ]; then
# rm -r "$DESTi"
# fi
#
# mpirun -np $use_cores Estia_baseline_ana1.out \
# --dir="$DESTi" --ncount=$ncount --gravitation \
# omegaa=$omega operationmode=0 theta_resolution=0.04 over_illumination=0.000 \
# sample=$sample sample_length=$sample_length sample_height=$sample_height \
# lambda_start=$lambda_start lambda_end=$lambda_end enable_gravity=1 enable_chopper=0 enable_polarizer=0 enable_analyzer=1
#
# sample_height=0.01
# DESTi=$DEST/ana2_10x10_10
# if [ -e "$DESTi" ]; then
# rm -r "$DESTi"
# fi
#
# mpirun -np $use_cores Estia_baseline.out \
# --dir="$DESTi" --ncount=$ncount --gravitation \
# omegaa=$omega operationmode=0 theta_resolution=0.04 over_illumination=0.000 \
# sample=$sample sample_length=$sample_length sample_height=$sample_height \
# lambda_start=$lambda_start lambda_end=$lambda_end enable_gravity=1 enable_chopper=0 enable_polarizer=0 enable_analyzer=2
DESTi=$DEST/pol_ref_10x0p5
if [ -e "$DESTi" ]; then
rm -r "$DESTi"
fi
# sample_height=0.001
# DESTi=$DEST/ana2_10x1_10
# if [ -e "$DESTi" ]; then
# rm -r "$DESTi"
# fi
#
# mpirun -np $use_cores Estia_baseline_ana2.out \
# --dir="$DESTi" --ncount=$ncount --gravitation \
# omegaa=$omega operationmode=0 theta_resolution=0.04 over_illumination=0.000 \
# sample=$sample sample_length=$sample_length sample_height=$sample_height \
# lambda_start=$lambda_start lambda_end=$lambda_end enable_gravity=1 enable_chopper=0 enable_polarizer=0 enable_analyzer=2
mpirun -np $use_cores --hostfile hostnames Estia_baseline.out \
--dir="$DESTi" --format=NeXus --ncount=$ncount \
sample=$sample sample_length=$sample_length sample_height=$sample_height \
lambda_start=$lambda_start lambda_end=$lambda_end enable_gravity=0 enable_chopper=0 \
enable_polarizer=0 enable_analyzer=0
###################### Reference and Ni-layer measurement 1x1mm² sample ####################
###################### Polarizers 10x5mm² and 10x0.5mm² VS ####################
ncount=2e8
sample_length=0.005
sample_height=0.01
DESTi=$DEST/pol_10x5
if [ -e "$DESTi" ]; then
rm -r "$DESTi"
fi
mpirun -np $use_cores --hostfile hostnames Estia_baseline.out \
--dir="$DESTi" --format=NeXus --ncount=$ncount \
sample=$sample sample_length=$sample_length sample_height=$sample_height \
lambda_start=$lambda_start lambda_end=$lambda_end enable_gravity=0 enable_chopper=0 \
enable_polarizer=1 enable_analyzer=0
ncount=1e9
sample_length=0.0005
sample_height=0.01
DESTi=$DEST/pol_10x0p5
if [ -e "$DESTi" ]; then
rm -r "$DESTi"
fi
mpirun -np $use_cores --hostfile hostnames Estia_baseline.out \
--dir="$DESTi" --format=NeXus --ncount=$ncount \
sample=$sample sample_length=$sample_length sample_height=$sample_height \
lambda_start=$lambda_start lambda_end=$lambda_end enable_gravity=0 enable_chopper=0 \
enable_polarizer=1 enable_analyzer=0
+41 -35
View File
@@ -2,45 +2,51 @@
#-*- coding: utf8 -*-
import sys, os
import numpy as np
from subprocess import call
from numpy import *
from scipy.optimize import leastsq
FOLDER='../results2/'
PREFIX='selene_geo_'
SEED=3798
ITEMS=1000
CALL='mpirun -np 6 Selene_geometry.out --dir=%s --ncount=4e8 enable_gravity=0 sample_length=0.003 '
CALL+='tx_1=%f tz_1=%f ry_1=%f rz_1=%f tx_2=%f tz_2=%f ry_2=%f rz_2=%f '
CALL='/afs/psi.ch/project/sinq/sl6-64/bin/mpirun -np $SLURM_NPROCS Estia_baseline.out '
CALL=+'--dir=../results/selene_geo --ncount=1e9 '
CALL=+'omegaa=1.0 sample=4 sample_length=0.001 sample_height=0.01 '
CALL=+'lambda_start=2.5 lambda_end=15.0 enable_gravity=1 enable_chopper=1 '
CALL=+'lambda_min=3.75 selene1_foot1y=%.4f selene1_foot2y=%.4f '
TX=[-0.001, 0.001]
TZ=[-0.005, 0.005]
RY=[-0.001, 0.001]
RZ=[-2.5, 2.5]
def B(x,w):
# box function with full width w
return float32(abs(x)<=(w/2.))
def G(x,I0,x0,sigma):
# Gaussian with intensity I0, center x0 and standard deviation sigma
return I0*exp(-0.5*(x-x0)**2/sigma**2)
def Intensity(x, p):
I0, x0, sigma, w=p
return convolve(B(xc,w), G(xc, I0, x0, sigma), mode='same')
def Beam(x,I0,x0,sigma,w):
I=I0*where((x-x0)<(-w/2.), exp(-0.5*(x-x0+w/2.)**2/sigma**2),
where((x-x0)>(w/2.), exp(-0.5*(x-x0-w/2.)**2/sigma**2), 1.))
return I
def residuals(p, x, y):
I0, x0, sigma, w=p
return y-Beam(x, I0, x0, sigma, w)
def FWHM(pi):
rng=xc[where(Beam(xc, *pi)>=(pi[0]/2.))[0]]
return rng[-1]-rng[0]
x=linspace(-10.0005, 10.0005, 20001)
xc=(x[1:]+x[:-1])/2.
def analyze(fi):
sim=mr.McSim(fi)
det=sim['tof_sample']
ignore,y=det.project1d('x', bins=x/1000.)
p,res=leastsq(residuals, (y.max(), (xc*y).sum()/y.sum(), 0.01, 0.1), (xc, y))
return p,det
if __name__=='__main__':
call('mcstas -t -o Selene_geometry.c Selene_geometry.instr'.split())
call('mpicc -O2 -o Selene_geometry.out Selene_geometry.c -lm -DUSE_MPI'.split())
np.random.seed(SEED)
fh=open(os.path.join(FOLDER, PREFIX+'sims.dat'), 'w')
fh.write('# index tx_1 tz_1 ry_1 rz_1 tx_2 tz_2 ry_2 rz_2\n')
fh.write('# seed = %i\n'%SEED)
for i in range(ITEMS):
tx_1=np.random.random()*(TX[1]-TX[0])+TX[0]
tx_2=np.random.random()*(TX[1]-TX[0])+TX[0]
tz_1=np.random.random()*(TZ[1]-TZ[0])+TZ[0]
tz_2=np.random.random()*(TZ[1]-TZ[0])+TZ[0]
ry_1=np.random.random()*(RY[1]-RY[0])+RY[0]
ry_2=np.random.random()*(RY[1]-RY[0])+RY[0]
rz_1=np.random.random()*(RZ[1]-RZ[0])+RZ[0]
rz_2=np.random.random()*(RZ[1]-RZ[0])+RZ[0]
ln=' '.join(['%f'%val for val in [tx_1, tz_1, ry_1, rz_1, tx_2, tz_2, ry_2, rz_2]])
fh.write(ln+'\n')
print i, ':', ln
call((CALL%(os.path.join(FOLDER, PREFIX+'%05i/'%i), tx_1, tz_1, ry_1, rz_1, tx_2, tz_2, ry_2, rz_2)).split())
fh.close()
pass
+40
View File
@@ -0,0 +1,40 @@
#!/bin/bash
#SBATCH -J mcEstia_1
#SBATCH -N 2
#SBATCH --ntasks-per-node=24
#SBATCH --time=1-00:00:00
#SBATCH --mail-type=fail
#SBATCH --mail-user=artur.glavic@psi.ch
#SBATCH -o stdout.log
#SBATCH -e stderr.log
#SBATCH --partition=ll_long
echo "Starting at `date`"
echo "Running on hosts: $SLURM_NODELIST"
echo "Running on $SLURM_NNODES nodes."
echo "Running on $SLURM_NPROCS processors."
echo "Current working directory is `pwd`"
/usr/bin/modulecmd tcsh load mcstas
bash compile_if_needed.sh
for y1 in $(seq -0.009 0.001 0.009)
do
for y2 in $(seq -0.009 0.001 0.009) $(seq -0.09 0.01 -0.01) $(seq 0.01 0.01 0.09) $(seq -0.5 0.1 -0.1) $(seq 0.1 0.1 0.5)
do
printf -v Y1 %.3f $y1
printf -v Y2 %.3f $y2
echo $Y1 $Y2
/afs/psi.ch/project/sinq/sl6-64/bin/mpirun -np $SLURM_NPROCS Estia_baseline.out \
--dir=../results/selene1_geo_$Y1\_$Y2 --ncount=1e9 \
omegaa=1.0 sample=4 sample_length=0.00005 sample_height=0.01 \
lambda_start=2.5 lambda_end=15.0 enable_gravity=1 enable_chopper=1 \
lambda_min=3.75 selene1_foot1y=$Y1 selene1_foot2y=$Y2
done
done
echo "Program finished with exit code $? at: `date`"
+38
View File
@@ -0,0 +1,38 @@
#!/bin/tcsh
#SBATCH -J mcEstia_S
#SBATCH -N 2
#SBATCH --ntasks-per-node=24
#SBATCH --time=2:00:00
#SBATCH --mail-type=fail
#SBATCH --mail-user=artur.glavic@psi.ch
#SBATCH -o stdout.log
#SBATCH -e stderr.log
#SBATCH --partition=ll_short
echo "Starting at `date`"
echo "Running on hosts: $SLURM_NODELIST"
echo "Running on $SLURM_NNODES nodes."
echo "Running on $SLURM_NPROCS processors."
echo "Current working directory is `pwd`"
#module unload intel
#module load intel
module load mcstas
bash compile_if_needed.sh
mpirun -np $SLURM_NPROCS Estia_baseline.out --format=NeXus \
--dir=$*
# --dir=../results/$1 --ncount=1e9 \
# sample=1 enable_polarizer=1 \
# omegaa=6.0 sample_length=0.048 sample_height=0.01
echo "Program finished with exit code $? at: `date`"
echo "Compressing hdf5 file..."
h5repack -i ../results/$1/mccode.h5 -o ../results/$1/mccode.h5.c -f /entry1/data/detector_list_p_x_y_t_L_sx_sy_sz/events:GZIP=5 -l /entry1/data/detector_list_p_x_y_t_L_sx_sy_sz/events:CHUNK=3072x8
mv ../results/$1/mccode.h5.c ../results/$1/mccode.h5
echo "Scipt finished."
+40
View File
@@ -0,0 +1,40 @@
#!/bin/bash
#SBATCH -J mcEstia_2
#SBATCH -N 2
#SBATCH --ntasks-per-node=24
#SBATCH --time=1-00:00:00
#SBATCH --mail-type=fail
#SBATCH --mail-user=artur.glavic@psi.ch
#SBATCH -o stdout2.log
#SBATCH -e stderr2.log
#SBATCH --partition=ll_long
echo "Starting at `date`"
echo "Running on hosts: $SLURM_NODELIST"
echo "Running on $SLURM_NNODES nodes."
echo "Running on $SLURM_NPROCS processors."
echo "Current working directory is `pwd`"
/usr/bin/modulecmd tcsh load mcstas
bash compile_if_needed.sh
for y1 in $(seq -0.09 0.01 -0.01)
do
for y2 in $(seq -0.009 0.001 0.009) $(seq -0.09 0.01 -0.01) $(seq 0.01 0.01 0.09) $(seq -0.5 0.1 -0.1) $(seq 0.1 0.1 0.5)
do
printf -v Y1 %.3f $y1
printf -v Y2 %.3f $y2
echo $Y1 $Y2
/afs/psi.ch/project/sinq/sl6-64/bin/mpirun -np $SLURM_NPROCS Estia_baseline.out \
--dir=../results/selene1_geo_$Y1\_$Y2 --ncount=1e9 \
omegaa=1.0 sample=4 sample_length=0.00005 sample_height=0.01 \
lambda_start=2.5 lambda_end=15.0 enable_gravity=1 enable_chopper=1 \
lambda_min=3.75 selene1_foot1y=$Y1 selene1_foot2y=$Y2
done
done
echo "Program finished with exit code $? at: `date`"
+40
View File
@@ -0,0 +1,40 @@
#!/bin/bash
#SBATCH -J mcEstia_3
#SBATCH -N 2
#SBATCH --ntasks-per-node=24
#SBATCH --time=1-00:00:00
#SBATCH --mail-type=fail
#SBATCH --mail-user=artur.glavic@psi.ch
#SBATCH -o stdout3.log
#SBATCH -e stderr3.log
#SBATCH --partition=ll_long
echo "Starting at `date`"
echo "Running on hosts: $SLURM_NODELIST"
echo "Running on $SLURM_NNODES nodes."
echo "Running on $SLURM_NPROCS processors."
echo "Current working directory is `pwd`"
/usr/bin/modulecmd tcsh load mcstas
bash compile_if_needed.sh
for y1 in $(seq 0.01 0.01 0.09)
do
for y2 in $(seq -0.009 0.001 0.009) $(seq -0.09 0.01 -0.01) $(seq 0.01 0.01 0.09) $(seq -0.5 0.1 -0.1) $(seq 0.1 0.1 0.5)
do
printf -v Y1 %.3f $y1
printf -v Y2 %.3f $y2
echo $Y1 $Y2
/afs/psi.ch/project/sinq/sl6-64/bin/mpirun -np $SLURM_NPROCS Estia_baseline.out \
--dir=../results/selene1_geo_$Y1\_$Y2 --ncount=1e9 \
omegaa=1.0 sample=4 sample_length=0.00005 sample_height=0.01 \
lambda_start=2.5 lambda_end=15.0 enable_gravity=1 enable_chopper=1 \
lambda_min=3.75 selene1_foot1y=$Y1 selene1_foot2y=$Y2
done
done
echo "Program finished with exit code $? at: `date`"
+40
View File
@@ -0,0 +1,40 @@
#!/bin/bash
#SBATCH -J mcEstia_4
#SBATCH -N 2
#SBATCH --ntasks-per-node=24
#SBATCH --time=1-00:00:00
#SBATCH --mail-type=fail
#SBATCH --mail-user=artur.glavic@psi.ch
#SBATCH -o stdout4.log
#SBATCH -e stderr4.log
#SBATCH --partition=ll_long
echo "Starting at `date`"
echo "Running on hosts: $SLURM_NODELIST"
echo "Running on $SLURM_NNODES nodes."
echo "Running on $SLURM_NPROCS processors."
echo "Current working directory is `pwd`"
/usr/bin/modulecmd tcsh load mcstas
bash compile_if_needed.sh
for y1 in $(seq -0.5 0.1 -0.1)
do
for y2 in $(seq -0.009 0.001 0.009) $(seq -0.09 0.01 -0.01) $(seq 0.01 0.01 0.09) $(seq -0.5 0.1 -0.1) $(seq 0.1 0.1 0.5)
do
printf -v Y1 %.3f $y1
printf -v Y2 %.3f $y2
echo $Y1 $Y2
/afs/psi.ch/project/sinq/sl6-64/bin/mpirun -np $SLURM_NPROCS Estia_baseline.out \
--dir=../results/selene1_geo_$Y1\_$Y2 --ncount=1e9 \
omegaa=1.0 sample=4 sample_length=0.00005 sample_height=0.01 \
lambda_start=2.5 lambda_end=15.0 enable_gravity=1 enable_chopper=1 \
lambda_min=3.75 selene1_foot1y=$Y1 selene1_foot2y=$Y2
done
done
echo "Program finished with exit code $? at: `date`"
+40
View File
@@ -0,0 +1,40 @@
#!/bin/bash
#SBATCH -J mcEstia_5
#SBATCH -N 2
#SBATCH --ntasks-per-node=24
#SBATCH --time=1-00:00:00
#SBATCH --mail-type=fail
#SBATCH --mail-user=artur.glavic@psi.ch
#SBATCH -o stdout5.log
#SBATCH -e stderr5.log
#SBATCH --partition=ll_long
echo "Starting at `date`"
echo "Running on hosts: $SLURM_NODELIST"
echo "Running on $SLURM_NNODES nodes."
echo "Running on $SLURM_NPROCS processors."
echo "Current working directory is `pwd`"
/usr/bin/modulecmd tcsh load mcstas
bash compile_if_needed.sh
for y1 in $(seq 0.1 0.1 0.5)
do
for y2 in $(seq -0.009 0.001 0.009) $(seq -0.09 0.01 -0.01) $(seq 0.01 0.01 0.09) $(seq -0.5 0.1 -0.1) $(seq 0.1 0.1 0.5)
do
printf -v Y1 %.3f $y1
printf -v Y2 %.3f $y2
echo $Y1 $Y2
/afs/psi.ch/project/sinq/sl6-64/bin/mpirun -np $SLURM_NPROCS Estia_baseline.out \
--dir=../results/selene1_geo_$Y1\_$Y2 --ncount=1e9 \
omegaa=1.0 sample=4 sample_length=0.00005 sample_height=0.01 \
lambda_start=2.5 lambda_end=15.0 enable_gravity=1 enable_chopper=1 \
lambda_min=3.75 selene1_foot1y=$Y1 selene1_foot2y=$Y2
done
done
echo "Program finished with exit code $? at: `date`"