Added geometry correction. This refs #4303

mantidproject · Jan 25, 2012 · 68896da · 68896da
1 parent c68dd65
commit 68896da
Show file tree

Hide file tree

Showing 2 changed files with 70 additions and 37 deletions.
diff --git a/Code/Mantid/Framework/PythonAPI/PythonAlgorithms/RefLReduction.py b/Code/Mantid/Framework/PythonAPI/PythonAlgorithms/RefLReduction.py
@@ -167,7 +167,9 @@ def PyExec(self):
                                                     maxX=maxX,
                                                     maxY=maxY,
                                                     source_to_detector=dMD,
-                                                    theta=theta)
+                                                    sample_to_detector=dSD,
+                                                    theta=theta,
+                                                    geo_correction=False)
 
         #_tof_axis = mt2.readX(0)[:]
         ########## This was used to test the R(Q) 
@@ -234,14 +236,17 @@ def PyExec(self):
         #of interest along the x-axis (to avoid the frame effect)
         mt3_norm = wks_utility.createIntegratedWorkspace(mtd[ws_norm_histo_data], 
                                                          "IntegratedNormWks",
-                                            fromXpixel=Xrange[0],
-                                            toXpixel=Xrange[1],
-                                            fromYpixel=BackfromYpixel,
-                                            toYpixel=BacktoYpixel,
-                                            maxX=maxX,
-                                            maxY=maxY,
-                                            source_to_detector=dMD,
-                                            theta=theta)
+                                                         fromXpixel=Xrange[0],
+                                                         toXpixel=Xrange[1],
+                                                         fromYpixel=BackfromYpixel,
+                                                         toYpixel=BacktoYpixel,
+                                                         maxX=maxX,
+                                                         maxY=maxY,
+                                                         cpix=data_cpix,
+                                                         source_to_detector=dMD,
+                                                         sample_to_detector=dSD,
+                                                         theta=theta,
+                                                         geo_correction=False)
 
         Transpose(InputWorkspace='IntegratedNormWks',
                   OutputWorkspace='TransposedID')
@@ -273,6 +278,8 @@ def PyExec(self):
 
         # Normalization           
         SumSpectra(InputWorkspace="NormWks", OutputWorkspace="NormWks")
+
+
         #### divide data by normalize histo workspace
         Divide(LHSWorkspace='DataWks',
                RHSWorkspace='NormWks',

diff --git a/Code/Mantid/scripts/reduction/instruments/reflectometer/wks_utility.py b/Code/Mantid/scripts/reduction/instruments/reflectometer/wks_utility.py
@@ -97,48 +97,71 @@ def createIntegratedWorkspace(mt1, outputWorkspace,
                               fromXpixel, toXpixel,
                               fromYpixel, toYpixel,
                               maxX=304, maxY=256, 
-                              source_to_detector=None, theta=None):
+                              cpix=None,
+                              source_to_detector=None, 
+                              sample_to_detector=None,
+                              theta=None,
+                              geo_correction=False):
     """
         This creates the integrated workspace over the second pixel range (304 here) and
         returns the new workspace handle
     """
     _tof_axis = mt1.readX(0)[:]
 
-    if source_to_detector is not None and theta is not None:
-        _const = float(4) * math.pi * m * source_to_detector / h
-        _q_axis = 1e-10 * _const * math.sin(theta) / (_tof_axis*1e-6)
+    if geo_correction:
+
+        yrange = arange(toYpixel-fromYpixel+1) + fromYpixel
+        _q_axis = convertToRvsQWithCorrection(mt1, 
+                                              dMD=sample_to_detector,
+                                              theta=theta,
+                                              tof=_tof_axis, 
+                                              yrange=yrange, 
+                                              cpix=cpix)
+
+        print shape(_q_axis)
+
+    else:
 
-    _y_axis = zeros((maxY, len(_q_axis) - 1))
-    _y_error_axis = zeros((maxY, len(_q_axis) - 1))
+        if source_to_detector is not None and theta is not None:
+            _const = float(4) * math.pi * m * source_to_detector / h
+            _q_axis = 1e-10 * _const * math.sin(theta) / (_tof_axis*1e-6)
 
-    x_size = toXpixel - fromXpixel + 1 
-    x_range = arange(x_size) + fromXpixel
+        _y_axis = zeros((maxY, len(_q_axis) - 1))
+        _y_error_axis = zeros((maxY, len(_q_axis) - 1))
 
-    y_size = toYpixel - fromYpixel + 1
-    y_range = arange(y_size) + fromYpixel
+        x_size = toXpixel - fromXpixel + 1 
+        x_range = arange(x_size) + fromXpixel
 
-    for x in x_range:
-        for y in y_range:
-            _index = int((maxY) * x + y)
-            _y_axis[y, :] += mt1.readY(_index)[:]
-            _y_error_axis[y, :] += ((mt1.readE(_index)[:]) * (mt1.readE(_index)[:]))
+        y_size = toYpixel - fromYpixel + 1
+        y_range = arange(y_size) + fromYpixel
+
+        for x in x_range:
+            for y in y_range:
+                _index = int((maxY) * x + y)
+                _y_axis[y, :] += mt1.readY(_index)[:]
+                _y_error_axis[y, :] += ((mt1.readE(_index)[:]) * (mt1.readE(_index)[:]))
 
-    _y_axis = _y_axis.flatten()
-    _y_error_axis = numpy.sqrt(_y_error_axis)
-    #plot_y_error_axis = _y_error_axis #for output testing only    -> plt.imshow(plot_y_error_axis, aspect='auto', origin='lower')
-    _y_error_axis = _y_error_axis.flatten()
+        _y_axis = _y_axis.flatten()
+        _y_error_axis = numpy.sqrt(_y_error_axis)
+        #plot_y_error_axis = _y_error_axis #for output testing only    -> plt.imshow(plot_y_error_axis, aspect='auto', origin='lower')
+        _y_error_axis = _y_error_axis.flatten()
 
-    #normalization by proton charge
-#    _y_axis /= (proton_charge * 1e-12)
+        #normalization by proton charge
+        #_y_axis /= (proton_charge * 1e-12)
 
-    _q_axis = _q_axis[::-1]
-    _y_axis = _y_axis[::-1]
-    _y_error_axis = _y_error_axis[::-1]
+        _q_axis = _q_axis[::-1]
+        _y_axis = _y_axis[::-1]
+        _y_error_axis = _y_error_axis[::-1]
 
-    CreateWorkspace(OutputWorkspace=outputWorkspace, 
-                    DataX=_q_axis, DataY=_y_axis, DataE=_y_error_axis, Nspec=maxY,
-                    UnitX="MomentumTransfer",ParentWorkspace=mt1)
-    mt2 = mtd[outputWorkspace]
+        CreateWorkspace(OutputWorkspace=outputWorkspace, 
+                        DataX=_q_axis, DataY=_y_axis, DataE=_y_error_axis, Nspec=maxY,
+                        UnitX="MomentumTransfer",ParentWorkspace=mt1)
+        mt2 = mtd[outputWorkspace]
+
+
+
+
+
     return mt2
 
 def angleUnitConversion(value, from_units='degree', to_units='rad'):
@@ -269,6 +292,9 @@ def convertToRvsQWithCorrection(mt, dMD=-1,theta=-1,tof=None, yrange=None, cpix=
             theta: angle of detector
     """
 
+    h = 6.626e-34 #m^2 kg s^-1
+    m = 1.675e-27 #kg
+
     sample = mt.getInstrument().getSample()
     source = mt.getInstrument().getSource()
     dSM = sample.getDistance(source)