Changeset 580ee7a in sasview

Timestamp:

Aug 5, 2011 10:18:09 AM (13 years ago)

Author:

Jae Cho <jhjcho@…>

Branches:

master, ESS_GUI, ESS_GUI_Docs, ESS_GUI_batch_fitting, ESS_GUI_bumps_abstraction, ESS_GUI_iss1116, ESS_GUI_iss879, ESS_GUI_iss959, ESS_GUI_opencl, ESS_GUI_ordering, ESS_GUI_sync_sascalc, costrafo411, magnetic_scatt, release-4.1.1, release-4.1.2, release-4.2.2, release_4.0.1, ticket-1009, ticket-1094-headless, ticket-1242-2d-resolution, ticket-1243, ticket-1249, ticket885, unittest-saveload

Children:

148f0dd

Parents:

1bf1bc0

Message:

extended rescal to tof

Location:

calculator

Files:

• instrument.py (modified) (8 diffs)
• resolution_calculator.py (modified) (23 diffs)
• test/utest_resolution_calculator.py (modified) (1 diff)

calculator/instrument.py

@@ -17 +17 @@
 
 _MIN = 0.0
-_MAX = 30.0
+_MAX = 50.0
 _INTENSITY = 368428
 _WAVE_LENGTH = 6.0
 _WAVE_SPREAD = 0.125
 _MASS = 1.67492729E-24 #[gr]
-X_VAL = numpy.linspace(0, 30, 310)
-_LAMBDA_ARRAY = [[0, 5, 10, 15, 20, 25, 30], [1, 1, 1, 1, 1, 1, 1]]
+_LAMBDA_ARRAY = [[0, 1e+16],[_INTENSITY, _INTENSITY]]
 
 class Aperture(object):
@@ -152 +151 @@
             self.distance = distance
             validate(distance[0])
-    
-        
+
+            
 class Neutron(object):
     """
@@ -160 +159 @@
     def __init__(self):
         
-        # intensity in counts/sec
-        self.intensity = _INTENSITY
         # neutron mass in cgs unit
         self.mass = _MASS
-        # wavelength spectrum 
-        self.spectrum = []
+
         # wavelength
         self.wavelength = _WAVE_LENGTH
         # wavelength spread (FWHM)
         self.wavelength_spread = _WAVE_SPREAD
-        # mean wavelength
-        self.mean = _WAVE_LENGTH
-        # wavelength = distribution after velocity selector
-        self.peak = []
-        # std of of the spectrum
-        self.std = None
+        # wavelength spectrum 
+        self.spectrum = self.get_default_spectrum()        
+        # intensity in counts/sec
+        self.intensity = numpy.interp(self.wavelength, 
+                                      self.spectrum[0],
+                                      self.spectrum[1],
+                                      0.0,
+                                      0.0)
         # min max range of the spectrum
-        self.min = _MIN
-        self.max = _MAX
-        # x-range of the spectrum
-        self.x_val = X_VAL
-        # default distribution function
-        # ex., 'gaussian', or an array
-        self.func = _LAMBDA_ARRAY
+        self.min = min(self.spectrum[0])
+        self.max = max(self.spectrum[0])
+        # wavelength band
+        self.band = [self.min, self.max]
+
         # default unit of the thickness
         self.wavelength_unit = 'A'
-        
+    
+    def set_full_band(self):
+        """
+        set band to default value
+        """
+        self.band = self.spectrum
+    def set_spectrum(self, spectrum):  
+        """
+        Set spectrum
+        
+        :param spectrum: numpy array
+        """
+        self.spectrum = spectrum 
+        self.setup_spectrum()
+         
+    def setup_spectrum(self):
+        """
+        To set the wavelength spectrum, and intensity, assumes 
+        wavelength is already within the spectrum
+        """
+        spectrum = self.spectrum
+        intensity = numpy.interp(self.wavelength, 
+                                      spectrum[0],
+                                      spectrum[1],
+                                      0.0,
+                                      0.0)
+        self.set_intensity(intensity)
+        # min max range of the spectrum
+        self.min = min(self.spectrum[0])
+        self.max = max(self.spectrum[0])
+        # set default band
+        self.set_band([self.min,self.max])
+        
+    def set_band(self, band=[]):
+        """
+        To set the wavelength band
+        
+        :param band: array of [min, max]
+        """
+        # check if the wavelength is in range
+        if min(band) < self.min or\
+                max(band) > self.max:
+            raise 
+        self.band = band
+
+         
     def set_intensity(self, intensity = 368428):
         """
@@ -198 +239 @@
         Sets the wavelength
         """
+        # check if the wavelength is in range
+        if wavelength < min(self.band) or\
+                wavelength > max(self.band):
+            raise 
         self.wavelength = wavelength
         validate(wavelength)
+        self.intensity = numpy.interp(self.wavelength, 
+                                  self.spectrum[0],
+                                  self.spectrum[1],
+                                  0.0,
+                                  0.0)
+
 
     def set_mass(self, mass = _MASS):
@@ -245 +296 @@
         """
         return self.wavelength
-    
-    def _set_mean(self):
-        """
-        To get mean value of wavelength
-        """
-        mean_value = numpy.mean(self.peak[0]*
-                                self.peak[1])
-        self.mean = mean_value
-
-        
-    def get_mean_peakvalue(self):
-        """
-        To get mean value of wavelength
-        """
-        mean_value = numpy.mean(self.peak[1])
-        return mean_value
             
     def get_spectrum(self):
@@ -267 +302 @@
         """
         return self.spectrum
+    def get_default_spectrum(self):
+        """
+        get default spectrum
+        """
+        return numpy.array(_LAMBDA_ARRAY)
+    def get_band(self):
+        """
+        To get the wavelength band
+        """
+        return self.band 
     
     def plot_spectrum(self):
@@ -275 +320 @@
         try:
             import matplotlib.pyplot as plt
-            plt.plot(self.x_val, self.spectrum, linewidth = 2, color = 'r')
+            plt.plot(self.spectrum[0], self.spectrum[1], linewidth = 2, color = 'r')
             plt.legend(['Spectrum'], loc = 'best')
             plt.show()
@@ -281 +326 @@
             raise RuntimeError, "Can't import matplotlib required to plot..."
         
-    def plot_peak(self):
-        """
-        To plot the wavelength peak spactrum
-        : requirment: matplotlib.pyplot
-        """
-        try:
-            min = self.mean * (1 - self.wavelength_spread)
-            max = self.mean * (1 + self.wavelength_spread)
-            x_val =   numpy.linspace(min, max, 310)
-            import matplotlib.pyplot as plt
-            plt.plot(x_val, self.peak, linewidth = 2, color = 'r')
-            plt.legend(['peak'], loc = 'best')
-            plt.show()
-        except:
-            raise RuntimeError, "Can't import matplotlib required to plot..."
-
+    
+class TOF(Neutron):
+    """
+    TOF: make list of wavelength and wave length spreads
+    """
+    def __init__(self):
+        """
+        Init
+        """
+        Neutron.__init__(self)
+        #self.switch = switch
+        self.wavelength_list = [self.wavelength]
+        self.wavelength_spread_list = [self.wavelength_spread]
+        self.intensity_list = self.get_intensity_list()
+    
+    def get_intensity_list(self):
+        """
+        get list of the intensity wrt wavelength_list
+        """
+        out = numpy.interp(self.wavelength_list, 
+                                      self.spectrum[0],
+                                      self.spectrum[1],
+                                      0.0,
+                                      0.0)
+        return out
+    
+    def get_wave_list(self):
+        """
+        Get wavelength and wavelength_spread list
+        """
+        return self.wavelength_list, self.wavelengthspread_list
+    
+    def set_wave_list(self, wavelength=[]):
+        """
+        Set wavelength list
+        
+        :param wavelength: list of wavelengths
+        """
+        self.wavelength_list = wavelength
+        
+    def set_wave_spread_list(self, wavelength_spread=[]):
+        """
+        Set wavelength_spread list
+        
+        :param wavelength_spread: list of wavelength spreads
+        """
+        self.wavelengthspread_list = wavelength_spread
+        
+        
 def validate(value = None):
     """

calculator/resolution_calculator.py

@@ -6 +6 @@
 from instrument import Sample
 from instrument import Detector
-from instrument import Neutron
+from instrument import TOF as Neutron
 from instrument import Aperture
 # import math stuffs
@@ -36 +36 @@
         # 2d image of the resolution
         self.image = []
+        self.image_lam = []
         # resolutions
         # lamda in r-direction
@@ -55 +56 @@
         self.detector_qy_min = -0.3
         self.detector_qy_max = 0.3
+        # possible max qrange
+        self.qxmin_limit = 0
+        self.qxmax_limit = 0
+        self.qymin_limit = 0
+        self.qymax_limit = 0
+        # set sigmas
+        self.sigma_1 = 0
+        self.sigma_lamd = 0
+        self.sigma_2 = 0
+        self.sigma_1d = 0
         # plots
         self.plot = None
@@ -70 +81 @@
         self.detector_size = []
         # get all the values of the instrumental parameters
+        #self.intensity = self.get_intensity()
+        #self.wavelength = self.get_wavelength()
+        #self.wavelength_spread = self.get_wavelength_spread()
         self.get_all_instrument_params()
+        # max q range for all lambdas
+        self.qxrange = []
+        self.qyrange = []
         
     def compute_and_plot(self, qx_value, qy_value, qx_min, qx_max, 
@@ -79 +96 @@
         : qy_value: y component of q
         """
-        # compute 2d resolution
-        _, _, sigma_1, sigma_2, sigma_r = \
-                            self.compute(qx_value, qy_value, coord)
-        # make image
-        image = self.get_image(qx_value, qy_value, sigma_1, sigma_2, sigma_r, 
-                               qx_min, qx_max, qy_min, qy_max, coord)
+        # make sure to update all the variables need.
+        # except lambda, dlambda, and intensity
+        self.get_all_instrument_params()
+        # wavelength etc.
+        lamda_list, dlamb_list = self.get_wave_list()
+        intens_list = []#self.get_intensity_list()
+
+        sig1_list = []
+        sig2_list = []
+        sigr_list = []
+        sigma1d_list = []
+        num_lamda = len(lamda_list)
+        for num in range(num_lamda):
+            lam = lamda_list[num]
+            # wavelength spread
+            dlam = dlamb_list[num]
+            intens = self.setup_tof(lam, dlam)
+            intens_list.append(intens)
+            # save q min max
+            #qx_min = 0
+            #qx_max = 0
+            #qy_min = 0
+            #qy_max = 0
+            # compute 2d resolution
+            _, _, sigma_1, sigma_2, sigma_r, sigma1d = \
+                            self.compute(lam, dlam, qx_value, qy_value, coord)
+            # make image
+            image = self.get_image(qx_value, qy_value, sigma_1, sigma_2, 
+                            sigma_r, qx_min, qx_max, qy_min, qy_max,
+                            coord, False)
+            
+            # Non tof mode to be speed up 
+            #if num_lamda < 2:
+            #    return self.plot_image(image)
+            
+            if qx_min > self.qx_min:
+                qx_min = self.qx_min
+            if qx_max < self.qx_max:
+                qx_max = self.qx_max
+            if qy_min > self.qy_min:
+                qy_min = self.qy_min
+            if qy_max < self.qy_max:
+                qy_max = self.qy_max
+                
+            # set max qranges
+            self.qxrange = [qx_min, qx_max]
+            self.qyrange = [qy_min, qy_max]
+            
+            sig1_list.append(sigma_1)
+            sig2_list.append(sigma_2)
+            sigr_list.append(sigma_r)
+            sigma1d_list.append(sigma1d)
+        # redraw image in global 2d q-space.   
+        self.image_lam = []
+        total_intensity = 0
+        sigma_1 = 0
+        sigma_r = 0
+        sigma_2 = 0
+        sigma1d = 0
+        for ind in range(num_lamda):
+            lam = lamda_list[ind]
+            dlam = dlamb_list[ind]
+            intens = self.setup_tof(lam, dlam)
+            out = self.get_image(qx_value, qy_value, sig1_list[ind], 
+                                   sig2_list[ind], sigr_list[ind], 
+                                   qx_min, qx_max, qy_min, qy_max, coord)
+            # this is the case of q being outside the detector
+            #if numpy.all(out==0.0):
+            #    continue
+            image = out
+            # set sigmas
+            sigma_1 += sig1_list[ind] * self.intensity
+            sigma_r += sigr_list[ind] * self.intensity
+            sigma_2 += sig2_list[ind] * self.intensity
+            sigma1d += sigma1d_list[ind] * self.intensity
+            total_intensity += self.intensity
+        if total_intensity != 0:
+            image_out = image / total_intensity
+            sigma_1 = sigma_1 / total_intensity
+            sigma_r = sigma_r / total_intensity
+            sigma_2 = sigma_2 / total_intensity
+            sigma1d = sigma1d / total_intensity
+            # set sigmas
+            self.sigma_1 = sigma_1
+            self.sigma_lamd = sigma_r
+            self.sigma_2 = sigma_2
+            self.sigma_1d = sigma1d
+            # rescale 
+            max_im_val = 1 #image_out.max()
+            if max_im_val > 0:
+                image_out /= max_im_val
+        else:
+            image_out = image * 0.0
+            # Don't calculate sigmas nor set self.sigmas!
+            sigma_1 = 0
+            sigma_r = 0
+            sigma_2 = 0
+            sigma1d = 0
+        if len(self.image) > 0:
+            self.image += image_out
+        else:
+            self.image = image_out
+        
         # plot image
-        return self.plot_image(image) 
-
-    def compute(self, qx_value, qy_value, coord = 'cartesian'):
+        return self.plot_image(self.image) 
+    
+    def setup_tof(self, wavelength, wavelength_spread):
+        """
+        Setup all parameters in instrument
+        
+        : param ind: index of lambda, etc
+        """
+
+        # set wave.wavelength
+        self.set_wavelength(wavelength)
+        self.set_wavelength_spread(wavelength_spread)
+        self.intensity = self.wave.get_intensity()
+        
+        if wavelength == 0:
+            msg = "Can't compute the resolution: the wavelength is zero..." 
+            raise RuntimeError, msg
+        return self.intensity
+        
+    def compute(self, wavelength, wavelength_spread, qx_value, qy_value, 
+                coord = 'cartesian'):
         """
         Compute the Q resoltuion in || and + direction of 2D
@@ -95 +227 @@
         """
         coord = 'cartesian'
-        # make sure to update all the variables need.
-        self.get_all_instrument_params()
-        # wavelength
-        lamb = self.wavelength
-   
-        if lamb == 0:
-            msg = "Can't compute the resolution: the wavelength is zero..." 
-            raise RuntimeError, msg
-        # wavelength spread
-        lamb_spread = self.wavelength_spread
+        lamb = wavelength
+        lamb_spread = wavelength_spread
+        
         # Find polar values
         qr_value, phi = self._get_polar_value(qx_value, qy_value)
@@ -136 +261 @@
         lp_cor = (l_ssa * l_two) / (l_one + l_two)
         # the radial distance to the pixel from the center of the detector
-        radius = math.tan(theta)*l_two
+        radius = math.tan(theta) * l_two
         #Lp = l_one*l_two/(l_one+l_two)
         # default polar coordinate
@@ -161 +286 @@
                                           phi, 'radial', 'on')
         # for 1d
-        variance_1d_1 = sigma_1/2 +sigma_wave_1
+        variance_1d_1 = sigma_1/2 + sigma_wave_1
         # normalize
-        variance_1d_1 = knot*knot*variance_1d_1/12
+        variance_1d_1 = knot * knot * variance_1d_1 / 12
         
         # for 2d
         #sigma_1 += sigma_wave_1
         # normalize
-        sigma_1 = knot*sqrt(sigma_1/12)
-        sigma_r = knot*sqrt(sigma_wave_1/12)
+        sigma_1 = knot*sqrt(sigma_1 / 12)
+        sigma_r = knot*sqrt(sigma_wave_1 / 12)
         # sigma in the phi/y direction
         # for source apperture
@@ -192 +317 @@
         variance_1d_2 = sigma_2/2 +sigma_wave_2
         # normalize
-        variance_1d_2 = knot*knot*variance_1d_2/12
+        variance_1d_2 = knot*knot*variance_1d_2 / 12
         
         # for 2d
@@ -198 +323 @@
         #sigma_2 += sigma_wave_2
         # normalize
-        sigma_2 =  knot*sqrt(sigma_2/12)
-
+        sigma_2 =  knot * sqrt(sigma_2 / 12)
+        sigma1d = sqrt(variance_1d_1 + variance_1d_2)
         # set sigmas
         self.sigma_1 = sigma_1
         self.sigma_lamd = sigma_r
         self.sigma_2 = sigma_2
-        
-        self.sigma_1d = sqrt(variance_1d_1 + variance_1d_2)
-        return qr_value, phi, sigma_1, sigma_2, sigma_r
+        self.sigma_1d = sigma1d
+        return qr_value, phi, sigma_1, sigma_2, sigma_r, sigma1d
+    
+    def _within_detector_range(self,qx_value, qy_value):
+        """
+        check if qvalues are within detector range
+        """
+        # detector range
+        detector_qx_min = self.detector_qx_min
+        detector_qx_max = self.detector_qx_max
+        detector_qy_min = self.detector_qy_min
+        detector_qy_max = self.detector_qy_max
+        if self.qxmin_limit > detector_qx_min:
+            self.qxmin_limit = detector_qx_min
+        if self.qxmax_limit < detector_qx_max:
+            self.qxmax_limit = detector_qx_max
+        if self.qymin_limit > detector_qy_min:
+            self.qymin_limit = detector_qy_min
+        if self.qymax_limit < detector_qy_max:
+            self.qymax_limit = detector_qy_max
+        if qx_value < detector_qx_min or qx_value > detector_qx_max:
+            return False
+        if qy_value < detector_qy_min or qy_value > detector_qy_max:
+            return False
+        return True
     
     def get_image(self, qx_value, qy_value, sigma_1, sigma_2, sigma_r,
-                  qx_min, qx_max, qy_min, qy_max, coord = 'cartesian'): 
+                  qx_min, qx_max, qy_min, qy_max, 
+                  coord = 'cartesian', full_cal=True): 
         """
         Get the resolution in polar coordinate ready to plot
@@ -220 +368 @@
         # Get  qx_max and qy_max...
         output = self._get_detector_qxqy_pixels()
-        # Set qx_value/qy_value min/max
-        #qx_min = self.qx_min
-        #qx_max = self.qx_max
-        #qy_min = self.qy_min
-        #qy_max = self.qy_max
+       
         qr_value, phi = self._get_polar_value(qx_value, qy_value)
 
@@ -241 +385 @@
             self.qy_max = qy_max
             #raise ValueError, msg
-
+        if not full_cal:
+            return None
+ 
         # Make an empty graph in the detector scale
         dx_size = (self.qx_max - self.qx_min) / (1000 - 1)
@@ -269 +415 @@
             image = self._gaussian2d(q_1, q_2, qc_1, qc_2, 
                                      sigma_1, sigma_2, sigma_r)
+        # out side of detector
+        if not self._within_detector_range(qx_value, qy_value):
+            image *= 0.0
+            self.intensity = 0.0
+            #return self.image
+
         # Add it if there are more than one inputs.
-        if len(self.image) > 0:
-            self.image += image
-        else:
-            self.image = image
-
-        return self.image
+        if len(self.image_lam) > 0:
+            self.image_lam += image * self.intensity
+        else:
+            self.image_lam = image * self.intensity
+        
+        return self.image_lam
     
     def plot_image(self, image):
@@ -291 +443 @@
         # Max value of the image
         max = numpy.max(image)
+        qx_min, qx_max, qy_min, qy_max = self.get_detector_qrange()
+
         # Image
         im = plt.imshow(image, 
-                extent = [self.qx_min, self.qx_max, self.qy_min, self.qy_max])
+                extent = [qx_min, qx_max, qy_min, qy_max])
 
         # bilinear interpolation to make it smoother
@@ -446 +600 @@
         return self.wave.wavelength
     
+    def get_spectrum(self):
+        """
+        Get spectrum
+        """
+        return self.wave.spectrum   
+    
+    def get_default_spectrum(self):
+        """
+        Get default_spectrum
+        """
+        return self.wave.get_default_spectrum()    
+    
+    def get_spectrum(self):
+        """
+        Get _spectrum
+        """
+        return self.wave.get_spectrum() 
+         
     def get_wavelength_spread(self):
         """
@@ -475 +647 @@
         """
         return self.detector.pix_size
-    
+       
     def get_detector_size(self):
         """
@@ -505 +677 @@
         """
         self.wave.set_intensity(intensity)
-   
+        
+    def set_wave(self, wavelength):
+        """
+        Set wavelength list or wavelength
+        """
+        if wavelength.__class__.__name__ == 'list':
+            self.wave.set_wave_list(wavelength)
+        elif wavelength.__class__.__name__ == 'float':
+            self.wave.set_wave_list([wavelength])
+            #self.set_wavelength(wavelength)
+        else:
+            raise
+    
+    def set_wave_spread(self, wavelength_spread):
+        """
+        Set wavelength spread  or wavelength spread
+        """
+        if wavelength_spread.__class__.__name__ == 'list':
+            self.wave.set_wave_spread_list(wavelength_spread)
+        elif wavelength_spread.__class__.__name__ == 'float':
+            self.wave.set_wave_spread_list([wavelength_spread])
+            #self.set_wavelength_spread(wavelength_spread)
+        else:
+            raise
+        
     def set_wavelength(self, wavelength):
         """
         Set wavelength
         """
+        self.wavelength = wavelength
         self.wave.set_wavelength(wavelength)
-    
+        
+    def set_spectrum(self, spectrum):
+        """
+        Set spectrum
+        """
+        self.spectrum = spectrum
+        self.wave.set_spectrum(spectrum)  
+          
     def set_wavelength_spread(self, wavelength_spread):
         """
         Set wavelength spread
         """
+        self.wavelength_spread = wavelength_spread
         self.wave.set_wavelength_spread(wavelength_spread)
-    
+        
+    def set_wave_list(self, wavelength_list, wavelengthspread_list):
+        """
+        Set wavelength and its spread list
+        """
+        self.wave.set_wave_list(wavelength_list, wavelengthspread_list)
+    
+    def get_wave_list(self):
+        """
+        Set wavelength spread
+        """
+        return self.wave.get_wave_list()
+    
+    def get_intensity_list(self):
+        """
+        Set wavelength spread
+        """
+        return self.wave.get_intensity_list()
+           
     def set_source_aperture_size(self, size):
         """
@@ -592 +815 @@
         Get all instrumental parameters
         """
-        self.intensity = self.get_intensity()
-        self.wavelength = self.get_wavelength()
-        self.wavelength_spread = self.get_wavelength_spread()
+        #self.intensity = self.get_intensity()
+        #self.wavelength = self.get_wavelength()
+        #self.wavelength_spread = self.get_wavelength_spread()
         self.mass = self.get_neutron_mass()
+        self.spectrum = self.get_spectrum()
         self.source_aperture_size = self.get_source_aperture_size()
         self.sample_aperture_size = self.get_sample_aperture_size()
@@ -603 +827 @@
         self.sample2sample_distance = self.get_sample2sample_distance()
         self.sample2detector_distance = self.get_sample2detector_distance()
-
-
-          
+    
+    def get_detector_qrange(self):
+        """
+        get max detector q ranges
+        
+        : return: qx_min, qx_max, qy_min, qy_max tuple
+        """
+        if len(self.qxrange) != 2 or len(self.qyrange) != 2:
+            return None
+        qx_min = self.qxrange[0]
+        qx_max = self.qxrange[1]
+        qy_min = self.qyrange[0]
+        qy_max = self.qyrange[1]
+        
+        return qx_min, qx_max, qy_min, qy_max
+    
     def _rotate_z(self, x_value, y_value, theta= 0.0):
         """
@@ -744 +981 @@
         
         # wavelength
-        wavelength = self.wavelength
+        wavelength = self.wave.wavelength
         # Gavity correction
         delta_y = self._get_beamcenter_drop() # in cm
@@ -839 +1076 @@
         self.qx_max += self._get_qx(0.5 * pix_x_size, 
                                     sample2detector_distance, wavelength)
+        
         # min and max values of detecter
         self.detector_qx_min = self.qx_min
@@ -844 +1082 @@
         self.detector_qy_min = self.qy_min
         self.detector_qy_max = self.qy_max
+        
+        
 
         # try to set it as a Data2D otherwise pass (not required for now)
@@ -916 +1156 @@
         unit_cm = 1e-08
         # Velocity of neutron in horizontal direction (~ actual velocity)
-        velocity = _PLANK_H / (self.mass * self.wavelength * unit_cm)
+        velocity = _PLANK_H / (self.mass * self.wave.wavelength * unit_cm)
         # Compute delta y
         delta_y = 0.5

calculator/test/utest_resolution_calculator.py

@@ -24 +24 @@
         self.cal.set_sample2detector_distance([1500])
         self.cal.set_wavelength_spread(0.1)
-        qr, phi, sigma_1, sigma_2, _ = self.cal.compute(0, 0, coord = 'polar')
+        # all instr. params for cal
+        self.cal.get_all_instrument_params()
+        qr, phi, sigma_1, sigma_2, _, sigma_1d = self.cal.compute(15, 0.1, 
+                                                        0, 0, coord = 'polar')
         sigma_1d = self.cal.sigma_1d