Source Code for Module functions.dchi2

  1  from Numeric import Float64, zeros 
  2   
  3 -class dchi2: 
  4 -        def __init__(self): 
  5                  "Function to create the chi-squared gradient." 
  6   
  7   
  8 -        def dchi2(self, params, diff_type, diff_params, model, relax_data, errors, print_flag=0): 
  9                  """Function to create the chi-squared gradient. 
 10   
 11                  Function arguments 
 12                  ~~~~~~~~~~~~~~~~~~ 
 13   
 14                  1:  params - a list containing the parameter values specific for the given model. 
 15                  The order of parameters must be as follows: 
 16                          m1 - {S2} 
 17                          m2 - {S2, te} 
 18                          m3 - {S2, Rex} 
 19                          m4 - {S2, te, Rex} 
 20                          m5 - {S2f, S2s, ts} 
 21                  2:  diff_type - string.  The diffusion tensor, ie 'iso', 'axial', 'aniso' 
 22                  3:  diff_params - array.  An array with the diffusion parameters 
 23                  4:  model - string.  The model 
 24                  5:  relax_data - array.  An array containing the experimental relaxation values. 
 25                  6:  errors - array.  An array containing the experimental errors. 
 26   
 27   
 28                  The chi-sqared gradient 
 29                  ~~~~~~~~~~~~~~~~~~~~~~~ 
 30   
 31                  Data structure:  self.data.dchi2 
 32                  Dimension:  1D, (parameters) 
 33                  Type:  Numeric array, Float64 
 34                  Dependencies:  self.data.ri, self.data.dri 
 35                  Required by:  None 
 36   
 37   
 38                  Formula 
 39                  ~~~~~~~ 
 40                                 _n_ 
 41                   dChi2         \   /  Ri - Ri()      dRi()  \  
 42                  -------  =  -2  >  | ----------  .  ------- | 
 43                  dthetaj        /__ \ sigma_i**2     dthetaj / 
 44                                 i=1 
 45   
 46                  where: 
 47                          Ri are the values of the measured relaxation data set. 
 48                          Ri() are the values of the back calculated relaxation data set. 
 49                          sigma_i are the values of the error set. 
 50   
 51                  """ 
 52   
 53                  # debug. 
 54                  if print_flag == 2: 
 55                          print "\n< dchi2 >" 
 56                          print "Params: " + `params` 
 57   
 58                  self.data.params = params 
 59                  self.data.diff_type = diff_type 
 60                  self.data.diff_params = diff_params 
 61                  self.data.model = model 
 62                  self.data.relax_data = relax_data 
 63                  self.data.errors = errors 
 64   
 65                  # Test to see if relaxation values and spectral density values have previously been calculated for the current parameter values, 
 66                  # If not, calculate the relaxation values which will then calculate the spectral density values, and store the current data and parameter 
 67                  # values in self.data.relax_test and self.data.gradient_test. 
 68                  test = [ self.data.relax_data[:], self.data.errors[:], self.data.params.tolist(), self.data.model ] 
 69                  if test != self.data.relax_test: 
 70                          self.init_param_types() 
 71                          self.Ri() 
 72                          self.data.relax_test = test[:] 
 73                  if test != self.data.gradient_test: 
 74                          self.init_param_types() 
 75                          self.dRi() 
 76                          self.data.gradient_test = test[:] 
 77   
 78                  # Initialise the chi-squared gradient. 
 79                  self.data.dchi2 = zeros((len(self.data.params)), Float64) 
 80   
 81                  # Calculate the chi-squared gradient. 
 82                  for i in range(len(self.data.relax_data)): 
 83                          if self.data.errors[i] != 0.0: 
 84                                  # Parameter independent terms. 
 85                                  a = -2.0 * (self.data.relax_data[i] - self.data.ri[i]) / (self.data.errors[i]**2) 
 86                                  for param in range(len(self.data.params)): 
 87                                          self.data.dchi2[param] = self.data.dchi2[param] + a * self.data.dri[param, i] 
 88                          else: 
 89                                  for param in range(len(self.params)): 
 90                                          self.data.dchi2[param] = 1e99 
 91                                  break 
 92   
 93                  # debug. 
 94                  if print_flag == 2: 
 95                          print "J(w):  " + `self.data.jw` 
 96                          print "dJ(w): " + `self.data.djw` 
 97                          print "Ri:    " + `self.data.ri` 
 98                          print "dRi:   " + `self.data.dri` 
 99                          print "dchi2: " + `self.data.dchi2` 
100                          print "\n" 
101   
102                  return self.data.dchi2 
103