Source Code for Module topo.base.arrayutil

  1  """ 
  2  General utility functions and classes for Topographica that require numpy. 
  3   
  4  $Id: arrayutil.py 11290 2010-07-27 13:42:46Z ceball $ 
  5  """ 
  6  __version__ = "$Revision: 11290 $" 
  7   
  8  import re 
  9   
 10  from numpy import sqrt,dot,arctan2,array2string,fmod,floor,array, \ 
 11       unravel_index,concatenate,set_printoptions,divide,maximum,minimum 
 12  from numpy import abs # pylint: disable-msg=W0622 
 13  from numpy import ufunc 
 14   
 15  import param 
 16   
 17  # Ask numpy to print even relatively large arrays by default 
 18  set_printoptions(threshold=200*200) 
 19   
 20   
 21 -def ufunc_script_repr(f,imports,prefix=None,settings=None): 
 22      """ 
 23      Return a runnable representation of the numpy ufunc f, and an 
 24      import statement for its module. 
 25      """ 
 26      # (could probably be generalized if required, because module is 
 27      # f.__class__.__module__) 
 28      imports.append('import numpy') 
 29      return 'numpy.'+f.__name__ 
 30   
 31  from param import parameterized 
 32  parameterized.script_repr_reg[ufunc]=ufunc_script_repr 
 33   
 34   
 35 -def L2norm(v): 
 36      """ 
 37      Return the L2 norm of the vector v. 
 38      """ 
 39      return sqrt(dot(v,v)) 
 40   
 41   
 42 -def norm(v,p=2): 
 43      """ 
 44      Returns the Lp norm of v, where p is arbitrary and defaults to 2. 
 45      """ 
 46      return (abs(v)**p).sum()**(1.0/p) 
 47   
 48   
 49 -def divisive_normalization(weights): 
 50      """Divisively normalize an array to sum to 1.0""" 
 51      s = weights.sum() 
 52      if s != 0: 
 53          factor = 1.0/s 
 54          weights *= factor 
 55   
 56   
 57 -def add_border(matrix,width=1,value=0.0): 
 58      """    
 59      Returns a new matrix consisting of the given matrix with a border 
 60      or margin of the given width filled with the given value. 
 61      """ 
 62      rows,cols = matrix.shape 
 63   
 64      hborder = array([ [value]*(cols+2*width) ]*width) 
 65      vborder = array([ [value]*width ] * rows) 
 66   
 67      temp = concatenate( (vborder,matrix,vborder), axis=1) 
 68      return concatenate( (hborder,temp,hborder) )      
 69   
 70   
 71 -def arg(z): 
 72      """ 
 73      Return the complex argument (phase) of z. 
 74      (z in radians.) 
 75      """ 
 76      z = z + complex(0,0)  # so that arg(z) also works for real z 
 77      return arctan2(z.imag, z.real) 
 78   
 79   
 80 -def octave_str(mat,name="mat",owner=""): 
 81      """ 
 82      Print the given Numpy matrix in Octave format, listing the given 
 83      matrix name and the object that owns it (if any). 
 84      """ 
 85      # This just prints the string version of the matrix and does search/replace 
 86      # to convert it; there may be a faster or easier way. 
 87      mstr=array2string(mat) 
 88      mstr=re.sub('\n','',mstr) 
 89      mstr=re.sub('[[]','',mstr) 
 90      mstr=re.sub('[]]','\n',mstr) 
 91      return ("# Created from %s %s\n# name: %s\n# type: matrix\n# rows: %s\n# columns: %s\n%s" % 
 92             (owner,name,name,mat.shape[0],mat.shape[1],mstr)) 
 93   
 94   
 95 -def octave_output(filename,mat,name="mat",owner=""): 
 96      """Writes the given matrix to a new file of the given name, in Octave format.""" 
 97      f = open(filename,'w') 
 98      f.write(octave_str(mat,name,owner)) 
 99      f.close() 
100   
101   
102 -def centroid(array_2D): 
103      """Return the centroid (center of gravity) for a 2D array.""" 
104   
105      rows,cols = array_2D.shape 
106      rsum=0 
107      csum=0 
108      rmass_sum=0 
109      cmass_sum=0 
110      for r in xrange(rows): 
111          row_sum = array_2D[r,:].sum() 
112          rsum += r*row_sum 
113          rmass_sum += row_sum 
114       
115      for c in xrange(cols): 
116          col_sum = array_2D[:,c].sum() 
117          csum += c*col_sum 
118          cmass_sum += col_sum 
119           
120      row_centroid= rsum/rmass_sum 
121      col_centroid= csum/cmass_sum 
122   
123      return row_centroid, col_centroid 
124   
125   
126 -def clip_lower(arr,lower_bound): 
127      """ 
128      In-place, one-sided version of numpy.clip(). 
129   
130      i.e. numpy.clip(arr,a_min=lower_bound,out=arr) if it existed. 
131      """ 
132      maximum(arr,lower_bound,arr) 
133   
134   
135 -def clip_upper(arr,upper_bound): 
136      """ 
137      In-place, one-sided version of numpy.clip(). 
138   
139      i.e. numpy.clip(arr,a_max=upper_bound,out=arr) if it existed. 
140      """ 
141      minimum(arr,upper_bound,arr) 
142   
143   
144 -def wrap(lower, upper, x): 
145      """ 
146      Circularly alias the numeric value x into the range [lower,upper). 
147   
148      Valid for cyclic quantities like orientations or hues. 
149      """ 
150      #I have no idea how I came up with this algorithm; it should be simplified. 
151      # 
152      # Note that Python's % operator works on floats and arrays; 
153      # usually one can simply use that instead.  E.g. to wrap array or 
154      # scalar x into 0,2*pi, just use "x % (2*pi)". 
155      range_=upper-lower 
156      return lower + fmod(x-lower + 2*range_*(1-floor(x/(2*range_))), range_) 
157   
158   
159 -def array_argmax(arr): 
160      "Returns the coordinates of the maximum element in the given array." 
161      return unravel_index(arr.argmax(),arr.shape) 
162   
163   
164   
165  # CB: Is this of general interest? Used in gcal.ty. 
166 -class DivideWithConstant(param.Parameterized): 
167      """ 
168      Divide two scalars or arrays with a constant (c) offset on the 
169      denominator to allow setting the gain or to avoid divide-by-zero 
170      issues.  The non-constant part of the denominator (y) is clipped 
171      to ensure that it has only positive values. 
172      """ 
173      c = param.Number(default=1.0) 
174       
175 -    def __call__(self, x, y): 
176          return divide(x,maximum(y,0)+self.c) 
177