WFM simulation analysis helpers

fix sbatch script for mcc
Add Estia Rapid Subdivision Addon (ERSA) WFM device
2018-09-16 15:46:16 +02:00 · 2018-09-16 14:36:08 +02:00 · 2018-09-16 14:12:01 +02:00
24 changed files with 1147 additions and 1412 deletions
@@ -8,7 +8,6 @@ old
 *.backup
 .project
 .settings
-.idea
 *.pyc
 *.pyo
 hosts
@@ -1,9 +0,0 @@
- #!/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
@@ -0,0 +1,93 @@
+#!/usr/bin/env python
+
+from numpy import *
+
+import mcstas_reader as mr
+
+l=arange(3,11,0.01)
+lp=(l[1:]+l[:-1])/2.
+
+t=arange(0.0, 0.18, 0.0001)
+tp=(t[1:]+t[:-1])/2.
+
+def get_pulses(fn, freq, offset=0.002):
+    det=mr.McSim(fn)['tof_detector']
+    pulses=[]
+    for i in range(int(freq/7)):
+        li, Ii=det.project1d('L', bins=l, newcols=[('tc', 't-L*0.00406-%f'%offset)], 
+                             fltr='(tc>%f)&(tc<%f)'%(i*1./freq, (i+1)*1./freq))
+        ti, Iti=det.project1d('t', bins=t, newcols=[('tc', 't-L*0.00406-%f'%offset)], 
+                              fltr='(tc>%f)&(tc<%f)'%(i*1./freq, (i+1)*1./freq))
+        pulses.append(((li[1:]+li[:-1])/2., Ii, (ti[1:]+ti[:-1])/2., Iti))
+    return det, pulses
+
+def get_filtered_pulses(det, pulses, freq, offset=0.0001, threashold=0.01):
+    Is=array([pi[3] for pi in pulses])
+    fpulses=[]
+    for i in range(1,int(freq/7)-1):
+        fltr=(Is[:i].sum(axis=0)+Is[i+1:].sum(axis=0))<(threashold*Is[i])
+        tf=tp[fltr]
+        if len(tf)<2:
+            continue
+        li, Ii=det.project1d('L', bins=l, fltr='(t>%f)&(t<%f)'%(tf[0], tf[-1]))
+        ti, Iti=det.project1d('t', bins=t, fltr='(t>%f)&(t<%f)'%(tf[0], tf[-1]))
+        fpulses.append(((li[1:]+li[:-1])/2., Ii, (ti[1:]+ti[:-1])/2., Iti))
+    return fpulses
+    
+
+def calc_usable(pulses, threashold=0.01):
+    Is=array([pi[3] for pi in pulses])
+    sI=0.
+    for i in range(1,Is.shape[0]):
+        sI+=Is[i][(Is[:i].sum(axis=0)+Is[i+1:].sum(axis=0))<(threashold*Is[i])].sum()
+    return sI
+
+def calc_sum(pulses):
+    Is=array([pi[3] for pi in pulses])
+    return Is.sum(axis=0)
+
+def tof_plot(det, pulses):
+    figure(figsize=(18,10))
+    for i, pi in enumerate(pulses):
+        semilogy(tp, pi[3], label='%i'%i)
+    ylim(100,1e5)
+
+def tof_th_plot(det, pulses, threashold=0.01):
+    figure(figsize=(18,10))
+    Is=array([pi[3] for pi in pulses])
+    for i, pi in enumerate(pulses):
+        fltr=(Is[:i].sum(axis=0)+Is[i+1:].sum(axis=0))<(threashold*Is[i])
+        semilogy(tp[fltr], pi[3][fltr], label='%i'%i)
+    ylim(100,1e5)
+
+def lambda_plot(det, pulses):
+    figure(figsize=(18,10))
+    for i, pi in enumerate(pulses):
+        semilogy(lp, pi[1])
+    semilogy(lp, det.project1d('L', bins=l)[1])
+    semilogy(lp, Ir)
+    ylim(100,1e6)
+    xlim(3,11)
+
+def export(fn, pulses):
+    It=calc_sum(pulses)
+    savetxt(fn+'_lpulses.dat', array([pules[0][0], It, Ir]+[pi[1] for pi in pulses]).T)
+    savetxt(fn+'_tpulses.dat', array([pules[0][2]]+[pi[3] for pi in pulses]).T)
+
+ref=mr.McSim('../results/ersa_ref/')['tof_detector']
+Ir=ref.project1d('L', bins=l)[1]
+
+freq_offsets={245: 0.00007, 280: 0.002, 308: 0.0015, 350: 0.003}
+rfreq_offsets={245: 0.0000, 280: 0.0002, 308: 0.0015, 350: 0.0007}
+
+results={}
+for key,val in freq_offsets.items():
+    results['%i'%key]=get_pulses('../results/ersa_%iHz/'%key, 2.*key, val)
+    try:
+        results['r300-%i'%key]=get_pulses('../results/ersa_r300_%iHz/'%key, 2.*key, val)
+    except:
+        continue
+
+for key,val in rfreq_offsets.items():
+        results['rev-%i'%key]=get_pulses('../results/ersa_reverse_%iHz/'%key, 2.*key, val)
+
@@ -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 list(self.info['data'].keys())
+    return self.info['data'].keys()

  def monitors(self):
    output=[]
-    for val in list(self.info['data'].values()):
+    for val in 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 list(self.info['data'].items()):
+    for key, val in 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 list(self.keys()):
+    elif item in 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 in ['Li', 'List'] and y_col=='p': # Detector_nD     
+    if x_col=='Li' 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 list(self.hdf['data'].keys()):
+    for item in self.hdf['data'].keys():
      node=self.hdf['data/'+item]
      info={}
-      for key, value in list(node.attrs.items()):
-        info[key.strip()]=value.strip().decode('ascii')
+      for key, value in node.attrs.items():
+        info[key.strip()]=value.strip()
      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 in ['Li', 'List'] and y_col=='p': # Detector_nD
+    if x_col=='Li' 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=list(map(float, self.info['xlimits'].split()))
+    limits=map(float, self.info['xlimits'].split())
    x=linspace(limits[0], limits[1], len(self.data))

    ax.errorbar(x, self.data, yerr=self.errors)
@@ -291,7 +291,7 @@ class Dataset(object):
      import pylab
      ax=pylab.gca()

-    limits=list(map(float, self.info['xylimits'].split()))
+    limits=map(float, self.info['xylimits'].split())

    if log:
      ax.imshow(self.data, origin='lower', extent=limits, aspect='auto', norm=LogNorm())
@@ -342,7 +342,7 @@ class TofData(Dataset):
    if fltr is None:
      I, x=histogram(columns[col], bins=bins, weights=w)
    else:
-      if isinstance(fltr, str):
+      if isinstance(fltr, basestring):
        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, str):
+      if isinstance(fltr, basestring):
        fltr=eval(fltr, globals(), columns)
      I, y, x=histogram2d(columns[ycol][fltr], columns[xcol][fltr],
                        bins=bins, weights=columns['p'][fltr])
@@ -52,9 +52,9 @@ DEFINE INSTRUMENT ESS_reflectometer_Estia
 (double omegaa = 1.2,  int sample = 0,  double sample_length   = 0.01, double sample_height   = 0.01,
 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 pol1_start=0.3, double pol1_angle=1.66, double pol2_start=0.3, double pol2_angle=1.66,
+int enable_chopper = 0, int enable_gravity=0, int enable_windows=1, 
+double wfm_freq=0.0, double wfm_shift=0.0, double wfm_radius=0.15,
+int enable_polarizer = 0, int enable_analyzer = 0,
 double source_power = 2,
 double selene1_foot1y =0.0, double selene1_foot2y = 0.0
 )
@@ -105,22 +105,14 @@ double velocity_max           ;     // m/s  neutron velocity of lambda_min


 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
@@ -131,7 +123,7 @@ 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
@@ -144,8 +136,8 @@ Theta2_analyzer2=atan((dist_ana_vfocus+analyzer2_start+analyzer2_length)/dist_an
 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);
+printf(" Selene 1 rotation  = %.4f deg\n", selene1_rot);
+printf(" Selene 1 shift  = %.1f mm\n", selene1_shift*1e3);

 %}
 TRACE
@@ -211,7 +203,7 @@ 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=3, Lmin = lambda_start, Lmax = lambda_end,
+    target_index=4, Lmin = lambda_start, Lmax = lambda_end,
    n_pulses = 1+enable_chopper, acc_power=source_power)
    AT (0, 0, 0) RELATIVE origin 

@@ -243,15 +235,6 @@ COMPONENT chopper = DiskChopper(radius=chopper_diameter/2.0, yheight=0.02,
    AT (0, 0, chopper_pos-2*NBOA_c) RELATIVE arm_virtual_source_beam
    

-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(
@@ -300,8 +283,8 @@ COMPONENT ac_slit = Slit(
 * Sample area *
 ***************/
    COMPONENT tof_sample = Monitor_nD(
-        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",
+        filename = "tof_sample",
+        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 list all",
        xwidth=0.05, yheight = 0.05)
        WHEN sample==4
        AT (0, 0, 0) RELATIVE arm_sample_beam
@@ -349,67 +332,36 @@ COMPONENT si_sample = Mirror(
    ROTATED (0, 90, 0) RELATIVE arm_sample


-
-
-
-/* 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 (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;
-    %}
+    AT (0, 0, 0) RELATIVE arm_detector
+    ROTATED (-selene_theta+(Theta1_analyzer1-Theta2_analyzer1)/2.0, 0, 0) RELATIVE arm_detector

 COMPONENT arm_analyzer2 = Arm()
-    AT (0, 0, 0) RELATIVE arm_detector2
-    ROTATED (-selene_theta+(Theta1_analyzer2-Theta2_analyzer2)/2.0, 0, 0) RELATIVE arm_detector2
+    AT (0, 0, 0) RELATIVE arm_detector
+    ROTATED (-selene_theta+(Theta1_analyzer2-Theta2_analyzer2)/2.0, 0, 0) RELATIVE arm_detector

 /* polarization analyser */
-COMPONENT analyzer1 = Polariser(enable_ref=1, d_substrate = 5e-4, reflect_d=0, reflect_u=0, lin=analyzer1_start, length=analyzer1_length,
+COMPONENT analyzer1 = Polariser(nIncRefr=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(enable_ref=1, d_substrate = 5e-4, reflect_d=0, reflect_u=0, lin=analyzer2_start, length=analyzer2_length,
+COMPONENT analyzer2 = Polariser(nIncRefr=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_analyzer2
-        ROTATED (0,0,0.0) RELATIVE arm_analyzer2
+        AT (0, 0.0, analyzer2_start) RELATIVE arm_analyzer
+        ROTATED (0,0,0.0) RELATIVE arm_analyzer

 /* detector */
    COMPONENT tof_detector = Monitor_nD(
-        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",
+        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",
        xwidth = 0.5, yheight = 1.0)
        WHEN sample!=4
        AT (0, 0, detector_arm+0.00001) RELATIVE arm_detector
-        ROTATED (0, 0, 0) RELATIVE arm_detector

 /***********************************************************************/

@@ -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 = Absorber(xmin=NBOA_side_x, xmax=NBOA_side_x+0.012,
+COMPONENT NBOA_side = AbsorberFixed(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 = Absorber(xmin=-0.1, xmax=0.0,
+COMPONENT E02_01_011_wedge = AbsorberFixed(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 = Absorber(xmin=-0.1, xmax=0.0,
        ALLOW_BACKPROP;    
        PROP_Z0;
    %}
-COMPONENT E02_01_012_wedge = Absorber(xmin=-0.1, xmax=0.0,
+COMPONENT E02_01_012_wedge = AbsorberFixed(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 = Absorber(xmin=-0.1, xmax=0.0,
+COMPONENT E02_01_021_wedge = AbsorberFixed(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 = Absorber(xmin=-0.1, xmax=0.0,
        ALLOW_BACKPROP;    
        PROP_Z0;
    %}
-COMPONENT E02_01_022_wedge = Absorber(xmin=-0.1, xmax=0.0,
+COMPONENT E02_01_022_wedge = AbsorberFixed(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 = Absorber(xmin=-0.1, xmax=0.0,
+COMPONENT E02_01_031_wedge = AbsorberFixed(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 = Absorber(xmin=-0.1, xmax=0.0,
        ALLOW_BACKPROP;    
        PROP_Z0;
    %}
-COMPONENT E02_01_032_wedge = Absorber(xmin=-0.1, xmax=0.0,
+COMPONENT E02_01_032_wedge = AbsorberFixed(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 = Absorber(xmin=-0.1, xmax=0.0,
        ALLOW_BACKPROP;    
        PROP_Z0;
    %}
-COMPONENT E02_01_033_wedge = Absorber(xmin=-0.1, xmax=0.0,
+COMPONENT E02_01_033_wedge = AbsorberFixed(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 = Absorber(xmin=NFGA_side_x1, xmax=NFGA_side_x1+0.008,
+COMPONENT NFGA_side1 = AbsorberFixed(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 = Absorber(xmin=NFGA_side_x1, xmax=NFGA_side_x1+0.008,
        ALLOW_BACKPROP;    
        PROP_Z0;
    %}
-COMPONENT NFGA_side2 = Absorber(xmin=NFGA_side_x2, xmax=NFGA_side_x2+0.008,
+COMPONENT NFGA_side2 = AbsorberFixed(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 = Absorber(xmin=NFGA_side_x2, xmax=NFGA_side_x2+0.008,
    %}
        
 // in-bunker feeder segments	
-COMPONENT E02_02_011_wedge = Absorber(xmin=-0.1, xmax=0.0,
+COMPONENT E02_02_011_wedge = AbsorberFixed(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 = Absorber(xmin=-0.1, xmax=0.0,
        ALLOW_BACKPROP;    
        PROP_Z0;
    %}
-COMPONENT E02_02_012_wedge = Absorber(xmin=-0.1, xmax=0.0,
+COMPONENT E02_02_012_wedge = AbsorberFixed(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 = Absorber(xmin=-0.1, xmax=0.0,
+COMPONENT E02_02_021_wedge = AbsorberFixed(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 = Absorber(xmin=-0.1, xmax=0.0,
        ALLOW_BACKPROP;    
        PROP_Z0;
    %}
-COMPONENT E02_02_022_wedge = Absorber(xmin=-0.1, xmax=0.0,
+COMPONENT E02_02_022_wedge = AbsorberFixed(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 = Absorber(xmin=-0.1, xmax=0.0,
+COMPONENT E02_02_031_wedge = AbsorberFixed(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 = Absorber(xmin=-0.1, xmax=0.0,
        ALLOW_BACKPROP;    
        PROP_Z0;
    %}
-COMPONENT E02_02_032_wedge = Absorber(xmin=-0.1, xmax=0.0,
+COMPONENT E02_02_032_wedge = AbsorberFixed(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 = Absorber(xmin=-0.1, xmax=0.0,
        ALLOW_BACKPROP;    
        PROP_Z0;
    %}
-COMPONENT E02_02_033_wedge = Absorber(xmin=-0.1, xmax=0.0,
+COMPONENT E02_02_033_wedge = AbsorberFixed(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
@@ -28,80 +28,20 @@ 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;
+    
+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
 %{
+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

-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);
+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
@@ -123,118 +63,106 @@ REMOVABLE COMPONENT source = Moderator(radius = 0.001, focus_xw = selene_entry+0
    AT (0,0,0) RELATIVE ISCS
    ROTATED (-1.25, 1.25, 0) RELATIVE ISCS

+COMPONENT arm_ersa = Arm()
+    AT (0, 0, 0) RELATIVE arm_selene2
+    ROTATED (0, -selene_theta, 0) RELATIVE arm_selene2
+

 /**************************************
 * Middle focus between Selene guides *
 **************************************/
-COMPONENT arm_polref = Arm()
+COMPONENT arm_polarizer = 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 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 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 ersa_01 = DiskChopper(radius=wfm_radius, yheight=0.05,
+    theta_0=120, phase=60+wfm_shift*180,
+    nu=wfm_freq,  nslit=2)
+    WHEN wfm_freq>0
+    AT (0, 0, -0.08) RELATIVE arm_ersa

-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 ersa_02 = DiskChopper(radius=wfm_radius, yheight=0.05,
+    theta_0=120, phase=75+wfm_shift*180,
+    nu=wfm_freq,  nslit=2)
+    WHEN wfm_freq>0
+    AT (0, 0, -0.07) RELATIVE arm_ersa

-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 ersa_03 = DiskChopper(radius=wfm_radius, yheight=0.05,
+    theta_0=120, phase=90+wfm_shift*180,
+    nu=wfm_freq,  nslit=2)
+    WHEN wfm_freq>0
+    AT (0, 0, -0.06) RELATIVE arm_ersa

-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 ersa_04 = DiskChopper(radius=wfm_radius, yheight=0.05,
+    theta_0=120, phase=105+wfm_shift*180,
+    nu=wfm_freq,  nslit=2)
+    WHEN wfm_freq>0
+    AT (0, 0, -0.05) RELATIVE arm_ersa

-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 ersa_05 = DiskChopper(radius=wfm_radius, yheight=0.05,
+    theta_0=120, phase=120+wfm_shift*180,
+    nu=wfm_freq,  nslit=2)
+    WHEN wfm_freq>0
+    AT (0, 0, -0.04) RELATIVE arm_ersa

-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 ersa_05b = DiskChopper(radius=wfm_radius, yheight=0.05,
+    theta_0=45, phase=-22.5+wfm_shift*180,
+    nu=wfm_freq,  nslit=2)
+    WHEN wfm_freq<0
+    AT (0, 0, -0.04) RELATIVE arm_ersa
+
+COMPONENT ersa_11b = DiskChopper(radius=wfm_radius, yheight=0.05,
+    theta_0=45, phase=22.5+wfm_shift*180,
+    nu=wfm_freq,  nslit=2)
+    WHEN wfm_freq<0
+    AT (0, 0, 0.04) RELATIVE arm_ersa
+
+COMPONENT ersa_11 = DiskChopper(radius=wfm_radius, yheight=0.05,
+    theta_0=120, phase=120+wfm_shift*180,
+    nu=wfm_freq,  nslit=2)
+    WHEN wfm_freq>0
+    AT (0, 0, 0.04) RELATIVE arm_ersa
+
+COMPONENT ersa_12 = DiskChopper(radius=wfm_radius, yheight=0.05,
+    theta_0=120, phase=135+wfm_shift*180,
+    nu=wfm_freq,  nslit=2)
+    WHEN wfm_freq>0
+    AT (0, 0, 0.05) RELATIVE arm_ersa
+
+COMPONENT ersa_13 = DiskChopper(radius=wfm_radius, yheight=0.05,
+    theta_0=120, phase=150+wfm_shift*180,
+    nu=wfm_freq,  nslit=2)
+    WHEN wfm_freq>0
+    AT (0, 0, 0.06) RELATIVE arm_ersa
+
+COMPONENT ersa_14 = DiskChopper(radius=wfm_radius, yheight=0.05,
+    theta_0=120, phase=165+wfm_shift*180,
+    nu=wfm_freq,  nslit=2)
+    WHEN wfm_freq>0
+    AT (0, 0, 0.07) RELATIVE arm_ersa
+
+COMPONENT ersa_15 = DiskChopper(radius=wfm_radius, yheight=0.05,
+    theta_0=120, phase=180+wfm_shift*180,
+    nu=wfm_freq,  nslit=2)
+    WHEN wfm_freq>0
+    AT (0, 0, 0.08) RELATIVE arm_ersa
+
+/* 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
+        


-
-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
@@ -102,7 +102,7 @@ REMOVABLE COMPONENT source = Moderator(radius = 0.001, focus_xw = selene_entry+0
 /* 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)
+    ymin = 0, ymax = selene_entry+0.005)
    AT (0, 0, selene_distance-0.01) RELATIVE arm_selene1

 COMPONENT block_before_selene_guide_1 = Absorber(
@@ -118,7 +118,7 @@ COMPONENT E02_03_01 = Elliptic_guide_gravity(
    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,
+    mright=0, mleft=scoating_01, mtop=scoating_01, mbottom=0,
    enableGravity=enable_gravity) 
    AT (0, 0, selene_distance) RELATIVE arm_selene1

@@ -127,7 +127,7 @@ COMPONENT E02_03_02 = Elliptic_guide_gravity(
    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,
+    mright=0, mleft=scoating_02, mtop=scoating_02, mbottom=0,
    enableGravity=enable_gravity) 
    AT (0, 0, selene_distance+1*selene_segment) RELATIVE arm_selene1

@@ -136,7 +136,7 @@ COMPONENT E02_03_03 = Elliptic_guide_gravity(
    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,
+    mright=0, mleft=scoating_03, mtop=scoating_03, mbottom=0,
    enableGravity=enable_gravity) 
    AT (0, 0, selene_distance+2*selene_segment) RELATIVE arm_selene1

@@ -145,7 +145,7 @@ COMPONENT E02_03_04 = Elliptic_guide_gravity(
    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,
+    mright=0, mleft=scoating_04, mtop=scoating_04, mbottom=0,
    enableGravity=enable_gravity) 
    AT (0, 0, selene_distance+3*selene_segment) RELATIVE arm_selene1

@@ -154,7 +154,7 @@ COMPONENT E02_03_05 = Elliptic_guide_gravity(
    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,
+    mright=0, mleft=scoating_05, mtop=scoating_05, mbottom=0,
    enableGravity=enable_gravity) 
    AT (0, 0, selene_distance+4*selene_segment) RELATIVE arm_selene1

@@ -163,7 +163,7 @@ COMPONENT E02_03_06 = Elliptic_guide_gravity(
    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,
+    mright=0, mleft=scoating_06, mtop=scoating_06, mbottom=0,
    enableGravity=enable_gravity) 
    AT (0, 0, selene_distance+5*selene_segment) RELATIVE arm_selene1

@@ -172,7 +172,7 @@ COMPONENT E02_03_07 = Elliptic_guide_gravity(
    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,
+    mright=0, mleft=scoating_07, mtop=scoating_07, mbottom=0,
    enableGravity=enable_gravity) 
    AT (0, 0, selene_distance+6*selene_segment) RELATIVE arm_selene1

@@ -181,7 +181,7 @@ COMPONENT E02_03_08 = Elliptic_guide_gravity(
    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,
+    mright=0, mleft=scoating_08, mtop=scoating_08, mbottom=0,
    enableGravity=enable_gravity) 
    AT (0, 0, selene_distance+7*selene_segment) RELATIVE arm_selene1

@@ -190,7 +190,7 @@ COMPONENT E02_03_09 = Elliptic_guide_gravity(
    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,
+    mright=0, mleft=scoating_09, mtop=scoating_09, mbottom=0,
    enableGravity=enable_gravity) 
    AT (0, 0, selene_distance+8*selene_segment) RELATIVE arm_selene1

@@ -199,7 +199,7 @@ COMPONENT E02_03_10 = Elliptic_guide_gravity(
    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,
+    mright=0, mleft=scoating_10, mtop=scoating_10, mbottom=0,
    enableGravity=enable_gravity) 
    AT (0, 0, selene_distance+9*selene_segment) RELATIVE arm_selene1

@@ -208,7 +208,7 @@ COMPONENT E02_03_11 = Elliptic_guide_gravity(
    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,
+    mright=0, mleft=scoating_11, mtop=scoating_11, mbottom=0,
    enableGravity=enable_gravity) 
    AT (0, 0, selene_distance+10*selene_segment) RELATIVE arm_selene1

@@ -217,7 +217,7 @@ COMPONENT E02_03_12 = Elliptic_guide_gravity(
    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,
+    mright=0, mleft=scoating_12, mtop=scoating_12, mbottom=0,
    enableGravity=enable_gravity) 
    AT (0, 0, selene_distance+11*selene_segment) RELATIVE arm_selene1

@@ -226,7 +226,7 @@ COMPONENT E02_03_13 = Elliptic_guide_gravity(
    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,
+    mright=0, mleft=scoating_13, mtop=scoating_13, mbottom=0,
    enableGravity=enable_gravity) 
    AT (0, 0, selene_distance+12*selene_segment) RELATIVE arm_selene1

@@ -235,7 +235,7 @@ COMPONENT E02_03_14 = Elliptic_guide_gravity(
    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,
+    mright=0, mleft=scoating_14, mtop=scoating_14, mbottom=0,
    enableGravity=enable_gravity) 
    AT (0, 0, selene_distance+13*selene_segment) RELATIVE arm_selene1

@@ -244,14 +244,14 @@ COMPONENT E02_03_15 = Elliptic_guide_gravity(
    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,
+    mright=0, mleft=scoating_15, mtop=scoating_15, mbottom=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)
+    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(
@@ -42,7 +42,7 @@ TRACE
 /* 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)
+        ymin = -selene_entry-0.01, ymax = 0)
        AT (0, 0, selene_distance-0.002) RELATIVE arm_selene2

 COMPONENT block_before_selene_guide_2 = Absorber(
@@ -59,7 +59,7 @@ COMPONENT E02_04_01 = Elliptic_guide_gravity(
    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,
+    mright=scoating_15, mleft=0, mtop=0, mbottom=scoating_15,
    enableGravity=enable_gravity) 
    AT (0, 0, selene_distance) RELATIVE arm_selene2

@@ -68,7 +68,7 @@ COMPONENT E02_04_02 = Elliptic_guide_gravity(
    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,
+    mright=scoating_14, mleft=0, mtop=0, mbottom=scoating_14,
    enableGravity=enable_gravity) 
    AT (0, 0, selene_distance+1*selene_segment) RELATIVE arm_selene2

@@ -77,7 +77,7 @@ COMPONENT E02_04_03 = Elliptic_guide_gravity(
    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,
+    mright=scoating_13, mleft=0, mtop=0, mbottom=scoating_13,
    enableGravity=enable_gravity) 
    AT (0, 0, selene_distance+2*selene_segment) RELATIVE arm_selene2

@@ -86,7 +86,7 @@ COMPONENT E02_04_04 = Elliptic_guide_gravity(
    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,
+    mright=scoating_12, mleft=0, mtop=0, mbottom=scoating_12,
    enableGravity=enable_gravity) 
    AT (0, 0, selene_distance+3*selene_segment) RELATIVE arm_selene2

@@ -95,7 +95,7 @@ COMPONENT E02_04_05 = Elliptic_guide_gravity(
    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,
+    mright=scoating_11, mleft=0, mtop=0, mbottom=scoating_11,
    enableGravity=enable_gravity) 
    AT (0, 0, selene_distance+4*selene_segment) RELATIVE arm_selene2

@@ -104,7 +104,7 @@ COMPONENT E02_04_06 = Elliptic_guide_gravity(
    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,
+    mright=scoating_10, mleft=0, mtop=0, mbottom=scoating_10,
    enableGravity=enable_gravity) 
    AT (0, 0, selene_distance+5*selene_segment) RELATIVE arm_selene2

@@ -113,7 +113,7 @@ COMPONENT E02_04_07 = Elliptic_guide_gravity(
    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,
+    mright=scoating_09, mleft=0, mtop=0, mbottom=scoating_09,
    enableGravity=enable_gravity) 
    AT (0, 0, selene_distance+6*selene_segment) RELATIVE arm_selene2

@@ -122,7 +122,7 @@ COMPONENT E02_04_08 = Elliptic_guide_gravity(
    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,
+    mright=scoating_08, mleft=0, mtop=0, mbottom=scoating_08,
    enableGravity=enable_gravity) 
    AT (0, 0, selene_distance+7*selene_segment) RELATIVE arm_selene2

@@ -131,7 +131,7 @@ COMPONENT E02_04_09 = Elliptic_guide_gravity(
    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,
+    mright=scoating_07, mleft=0, mtop=0, mbottom=scoating_07,
    enableGravity=enable_gravity) 
    AT (0, 0, selene_distance+8*selene_segment) RELATIVE arm_selene2

@@ -140,7 +140,7 @@ COMPONENT E02_04_10 = Elliptic_guide_gravity(
    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,
+    mright=scoating_06, mleft=0, mtop=0, mbottom=scoating_06,
    enableGravity=enable_gravity) 
    AT (0, 0, selene_distance+9*selene_segment) RELATIVE arm_selene2

@@ -149,7 +149,7 @@ COMPONENT E02_04_11 = Elliptic_guide_gravity(
    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,
+    mright=scoating_05, mleft=0, mtop=0, mbottom=scoating_05,
    enableGravity=enable_gravity) 
    AT (0, 0, selene_distance+10*selene_segment) RELATIVE arm_selene2

@@ -158,7 +158,7 @@ COMPONENT E02_04_12 = Elliptic_guide_gravity(
    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,
+    mright=scoating_04, mleft=0, mtop=0, mbottom=scoating_04,
    enableGravity=enable_gravity) 
    AT (0, 0, selene_distance+11*selene_segment) RELATIVE arm_selene2

@@ -167,7 +167,7 @@ COMPONENT E02_04_13 = Elliptic_guide_gravity(
    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,
+    mright=scoating_03, mleft=0, mtop=0, mbottom=scoating_03,
    enableGravity=enable_gravity) 
    AT (0, 0, selene_distance+12*selene_segment) RELATIVE arm_selene2

@@ -176,7 +176,7 @@ COMPONENT E02_04_14 = Elliptic_guide_gravity(
    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,
+    mright=scoating_02, mleft=0, mtop=0, mbottom=scoating_02,
    enableGravity=enable_gravity) 
    AT (0, 0, selene_distance+13*selene_segment) RELATIVE arm_selene2

@@ -185,14 +185,14 @@ COMPONENT E02_04_15 = Elliptic_guide_gravity(
    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,
+    mright=scoating_01, mleft=0, mtop=0, mbottom=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)
+    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(
@@ -1,91 +0,0 @@
-/*******************************************************************************
-*
-* 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
@@ -1,108 +0,0 @@
-/*******************************************************************************
-*         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
-
-
@@ -1,58 +0,0 @@
-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		
@@ -1,58 +0,0 @@
-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
@@ -1,6 +1,4 @@
 #!/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 ] \
@@ -8,6 +6,6 @@ if [ Estia_baseline.instr -nt Estia_baseline.out ] || [ ! -f 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 --trace
-	mpicc -O2 -o Estia_baseline.out Estia_baseline.c -lm -DUSE_MPI -DUSE_NEXUS -lNeXus
+	mcstas -o Estia_baseline.c Estia_baseline.instr
+	mpicc -O3 -o Estia_baseline.out Estia_baseline.c -lm -DUSE_MPI
 fi
@@ -1,3 +0,0 @@
-mcc05:72
-mcc06:72
-
@@ -1,83 +0,0 @@
-#!/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')
-            
@@ -1,7 +0,0 @@
-!#/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
@@ -1,77 +1,213 @@
 #!/bin/bash

 DEST=../results
-ncount=1e8
-use_cores=16
+ncount=2e7
+use_cores=4
 sample=4
-sample_length=0.005
+omega=0.8
+sample_length=0.01
 sample_height=0.01
 lambda_start=3.0
 lambda_end=32.0

-bash compile_if_needed.sh
+###################### 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

-###################### Reference measurement 10x5mm² and 10x0.5mm² VS ####################

-ncount=1e8
-sample_length=0.005
+###################### Reference and Ni-layer measurement 10x10mm² sample ####################
+# ncount=1e10
+sample_length=0.01
 sample_height=0.01
+lambda_start=3.5
+lambda_end=30.0
+sample=1

-DESTi=$DEST/pol_ref_10x5
+# 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
 if [ -e "$DESTi" ]; then
 	rm -r "$DESTi"
 fi

-mpirun  -np $use_cores --hostfile hostnames Estia_baseline.out \
-        --dir="$DESTi" --format=NeXus --ncount=$ncount \
+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=0 enable_chopper=0 \
-		enable_polarizer=0 enable_analyzer=0
-       
-ncount=5e8
-sample_length=0.0005
-sample_height=0.01
+        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

-DESTi=$DEST/pol_ref_10x0p5
-if [ -e "$DESTi" ]; then
-	rm -r "$DESTi"
-fi
+# 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

-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
+# 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
       

-###################### Polarizers 10x5mm² and 10x0.5mm² VS ####################
+###################### Reference and Ni-layer measurement 1x1mm² sample ####################

-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
@@ -0,0 +1,177 @@
+#!/bin/bash
+
+DEST=../results
+ncount=2e10
+use_cores=$SLURM_NPROCS
+sample=1
+omega=5.0
+sample_length=0.01
+sample_height=0.01
+lambda_start=3.00
+lambda_end=15.0
+
+######################## 1% WFM for 30cm and 60cm device #######################
+
+DESTi=$DEST/ersa_245Hz
+if [ -e "$DESTi" ]; then
+	rm -r "$DESTi"
+fi
+
+mpirun  -np $use_cores Estia_baseline.out \
+        --dir="$DESTi" --format=NeXuS --ncount=$ncount --gravitation \
+        omegaa=$omega operationmode=0 theta_resolution=0.04 over_illumination=0.001 \
+        sample=$sample sample_length=$sample_length sample_height=$sample_height \
+        enable_windows=0 wfm_freq=245 wfm_shift=0.0 \
+        lambda_start=$lambda_start lambda_end=$lambda_end enable_gravity=1 enable_chopper=1
+
+
+DESTi=$DEST/ersa_280Hz
+if [ -e "$DESTi" ]; then
+	rm -r "$DESTi"
+fi
+
+mpirun  -np $use_cores Estia_baseline.out \
+        --dir="$DESTi" --format=NeXuS --ncount=$ncount --gravitation \
+        omegaa=$omega operationmode=0 theta_resolution=0.04 over_illumination=0.001 \
+        sample=$sample sample_length=$sample_length sample_height=$sample_height \
+        enable_windows=0 wfm_freq=280 wfm_shift=0.0 \
+        lambda_start=$lambda_start lambda_end=$lambda_end enable_gravity=1 enable_chopper=1
+
+
+DESTi=$DEST/ersa_308Hz
+if [ -e "$DESTi" ]; then
+	rm -r "$DESTi"
+fi
+
+mpirun  -np $use_cores Estia_baseline.out \
+        --dir="$DESTi" --format=NeXuS --ncount=$ncount --gravitation \
+        omegaa=$omega operationmode=0 theta_resolution=0.04 over_illumination=0.001 \
+        sample=$sample sample_length=$sample_length sample_height=$sample_height \
+        enable_windows=0 wfm_freq=308 wfm_shift=0.0 \
+        lambda_start=$lambda_start lambda_end=$lambda_end enable_gravity=1 enable_chopper=1
+
+
+DESTi=$DEST/ersa_350Hz
+if [ -e "$DESTi" ]; then
+	rm -r "$DESTi"
+fi
+
+mpirun  -np $use_cores Estia_baseline.out \
+        --dir="$DESTi" --format=NeXuS --ncount=$ncount --gravitation \
+        omegaa=$omega operationmode=0 theta_resolution=0.04 over_illumination=0.001 \
+        sample=$sample sample_length=$sample_length sample_height=$sample_height \
+        enable_windows=0 wfm_freq=350 wfm_shift=0.0 \
+        lambda_start=$lambda_start lambda_end=$lambda_end enable_gravity=1 enable_chopper=1
+
+
+
+DESTi=$DEST/ersa_r300_245Hz
+if [ -e "$DESTi" ]; then
+	rm -r "$DESTi"
+fi
+
+mpirun  -np $use_cores Estia_baseline.out \
+        --dir="$DESTi" --format=NeXuS --ncount=$ncount --gravitation \
+        omegaa=$omega operationmode=0 theta_resolution=0.04 over_illumination=0.001 \
+        sample=$sample sample_length=$sample_length sample_height=$sample_height \
+        enable_windows=0 wfm_freq=245 wfm_shift=0.0 wfm_radius=0.3 \
+        lambda_start=$lambda_start lambda_end=$lambda_end enable_gravity=1 enable_chopper=1
+
+
+DESTi=$DEST/ersa_r300_280Hz
+if [ -e "$DESTi" ]; then
+	rm -r "$DESTi"
+fi
+
+mpirun  -np $use_cores Estia_baseline.out \
+        --dir="$DESTi" --format=NeXuS --ncount=$ncount --gravitation \
+        omegaa=$omega operationmode=0 theta_resolution=0.04 over_illumination=0.001 \
+        sample=$sample sample_length=$sample_length sample_height=$sample_height \
+        enable_windows=0 wfm_freq=280 wfm_shift=0.0 wfm_radius=0.3 \
+        lambda_start=$lambda_start lambda_end=$lambda_end enable_gravity=1 enable_chopper=1
+
+
+
+DESTi=$DEST/ersa_r300_308Hz
+if [ -e "$DESTi" ]; then
+	rm -r "$DESTi"
+fi
+
+mpirun  -np $use_cores Estia_baseline.out \
+        --dir="$DESTi" --format=NeXuS --ncount=$ncount --gravitation \
+        omegaa=$omega operationmode=0 theta_resolution=0.04 over_illumination=0.001 \
+        sample=$sample sample_length=$sample_length sample_height=$sample_height \
+        enable_windows=0 wfm_freq=308 wfm_shift=0.0 wfm_radius=0.3 \
+        lambda_start=$lambda_start lambda_end=$lambda_end enable_gravity=1 enable_chopper=1
+
+
+
+DESTi=$DEST/ersa_r300_350Hz
+if [ -e "$DESTi" ]; then
+	rm -r "$DESTi"
+fi
+
+mpirun  -np $use_cores Estia_baseline.out \
+        --dir="$DESTi" --format=NeXuS --ncount=$ncount --gravitation \
+        omegaa=$omega operationmode=0 theta_resolution=0.04 over_illumination=0.001 \
+        sample=$sample sample_length=$sample_length sample_height=$sample_height \
+        enable_windows=0 wfm_freq=350 wfm_shift=0.0 wfm_radius=0.3 \
+        lambda_start=$lambda_start lambda_end=$lambda_end enable_gravity=1 enable_chopper=1
+
+###################### 0.5% WFM ########################
+
+ncount=4e10
+DESTi=$DEST/ersa_reverse_280Hz
+if [ -e "$DESTi" ]; then
+	rm -r "$DESTi"
+fi
+
+mpirun  -np $use_cores Estia_baseline.out \
+        --dir="$DESTi" --format=NeXuS --ncount=$ncount \
+        omegaa=$omega operationmode=0 theta_resolution=0.04 over_illumination=0.001 \
+        sample=$sample sample_length=$sample_length sample_height=$sample_height \
+        enable_windows=0 wfm_freq=-280 wfm_shift=0.0 \
+        lambda_start=$lambda_start lambda_end=$lambda_end enable_gravity=0 enable_chopper=1
+
+
+DESTi=$DEST/ersa_reverse_308Hz
+if [ -e "$DESTi" ]; then
+	rm -r "$DESTi"
+fi
+
+mpirun  -np $use_cores Estia_baseline.out \
+        --dir="$DESTi" --format=NeXuS --ncount=$ncount \
+        omegaa=$omega operationmode=0 theta_resolution=0.04 over_illumination=0.001 \
+        sample=$sample sample_length=$sample_length sample_height=$sample_height \
+        enable_windows=0 wfm_freq=-308 wfm_shift=0.0 \
+        lambda_start=$lambda_start lambda_end=$lambda_end enable_gravity=0 enable_chopper=1
+
+
+DESTi=$DEST/ersa_reverse_350Hz
+if [ -e "$DESTi" ]; then
+	rm -r "$DESTi"
+fi
+
+mpirun  -np $use_cores Estia_baseline.out \
+        --dir="$DESTi" --format=NeXuS --ncount=$ncount \
+        omegaa=$omega operationmode=0 theta_resolution=0.04 over_illumination=0.001 \
+        sample=$sample sample_length=$sample_length sample_height=$sample_height \
+        enable_windows=0 wfm_freq=-350 wfm_shift=0.0 \
+        lambda_start=$lambda_start lambda_end=$lambda_end enable_gravity=0 enable_chopper=1
+
+
+##################### Reference run ###############################
+ncount=1e10
+
+DESTi=$DEST/ersa_ref
+if [ -e "$DESTi" ]; then
+	rm -r "$DESTi"
+fi
+
+mpirun  -np $use_cores Estia_baseline.out \
+        --dir="$DESTi" --format=NeXuS --ncount=$ncount --gravitation \
+        omegaa=$omega operationmode=0 theta_resolution=0.04 over_illumination=0.001 \
+        sample=$sample sample_length=$sample_length sample_height=$sample_height \
+        enable_windows=0 wfm_freq=0 \
+        lambda_start=$lambda_start lambda_end=$lambda_end enable_gravity=1 enable_chopper=1
+
@@ -0,0 +1,24 @@
+#!/bin/tcsh
+#SBATCH -J wfm_McEstia
+#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`"
+
+
+module load mcstas
+
+bash compile_if_needed.sh
+bash run_wfm.sh
@@ -1,38 +0,0 @@
-#!/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."
Author	SHA1	Message	Date
glavic_a	c14f216216	WFM simulation analysis helpers	2018-09-16 15:46:16 +02:00
glavic_a	a73791a813	fix sbatch script for mcc	2018-09-16 14:36:08 +02:00
glavic_a	d92e071e10	Add Estia Rapid Subdivision Addon (ERSA) WFM device	2018-09-16 14:12:01 +02:00