API tools faq
paste
Advertisement
chironex

Current _dissc.pyx

chironex
Apr 4th, 2011
123
0
Never
Add comment
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
Python 8.88 KB | None | 0 0
raw download clone embed print report
  1. # encoding: utf-8
  2. # cython: profile=False
  3. # filename: _dissc.pyx
  4.  
  5. '''
  6. Steve: Microtonal spectum optimizer? Cython Version
  7.  
  8. @author: Jordan
  9. '''
  10.  
  11. import numpy as np
  12. cimport numpy as np
  13. cimport cython
  14. import math as m
  15. from libc.math cimport exp as cExp
  16.  
  17. # DTYPE = np.float
  18. ctypedef np.float_t DTYPE_t
  19.  
  20. np.seterr(divide='raise', over='raise', under='ignore', invalid='raise')
  21.    
  22. """
  23. I define a timbre as the following 2d numpy array:
  24. [[f0, a0], [f1, a1], [f2, a2]...] where f describes the frequency
  25. of the given partial and a is its amplitude from 0 to 1. Phase is ignored.
  26. """
  27.  
  28. # MUST BE MODIFIED TO PROVIDE UNSORTED OUTPUT
  29. #Calculates the dissonance of a timbre applied to a scale.
  30. cpdef dissTimbreScale(np.ndarray[DTYPE_t,ndim=2] t,np.ndarray[DTYPE_t,ndim=1] s):
  31.     cdef DTYPE_t currDiss
  32.     currDiss = 0.0
  33.     cdef unsigned int i
  34.     for i in range(s.size):
  35.         currDiss += diss1Timbre(np.append(t[:,0],s[i]*t[:,0]), np.append(t[:,1],t[:,1]))
  36.     return currDiss
  37.  
  38. cdef inline DTYPE_t float_min(DTYPE_t a, DTYPE_t b): return a if a <= b else b
  39.  
  40. cdef inline DTYPE_t float_abs(DTYPE_t a): return a if a >= 0 else -a
  41.  
  42. # MODIFIED TO ACCEPT 2 1D UNSORTED INPUTS
  43. # Calculates the inherent dissonance of one timbres of the above form
  44. cpdef DTYPE_t diss1Timbre(np.ndarray[DTYPE_t,ndim=1] f, np.ndarray[DTYPE_t,ndim=1] a):
  45.     cdef DTYPE_t runningDiss1
  46.     runningDiss1 = 0.0
  47.     cdef unsigned int length = f.size
  48.     cdef unsigned int i
  49.     cdef unsigned int j
  50.     cdef DTYPE_t f1
  51.     cdef DTYPE_t f2
  52.     cdef DTYPE_t a1
  53.     cdef DTYPE_t a2
  54.     cdef DTYPE_t s
  55.     for i in range(length):
  56.         for j in range(i+1, length):
  57.            
  58.             if f[j] >= f[i]:
  59.                 f1 = f[i]
  60.                 f2 = f[j]
  61.                 a1 = float_abs(a[i])
  62.                 a2 = float_abs(a[j])
  63.             else:
  64.                 f2 = f[i]
  65.                 f1 = f[j]
  66.                 a2 = float_abs(a[i])
  67.                 a1 = float_abs(a[j])
  68.             s = 0.24/(0.021*f1 + 19)
  69.             runningDiss1 += (float_min(a1,a2) * (cExp(-3.5*s*(f2-f1)) - cExp(-5.75*s*(f2-f1)) ) )
  70.     return runningDiss1
  71.  
  72.  
  73.  
  74. # MODIFIED TO USE pDTS(freqs, amps, scale, i)
  75. # Calculates the gradient of the dissonance curve (does not consider amplitudes)
  76. cpdef np.ndarray[DTYPE_t, ndim=1] gradTimbreScale(np.ndarray[DTYPE_t, ndim=2] t,np.ndarray[DTYPE_t,ndim=1] s):
  77.     cdef np.ndarray[DTYPE_t, ndim=1] df
  78.     df = np.array([])
  79.     cdef unsigned int i
  80.     for i in range(1, t.shape[0]):
  81.         df = np.append(df,partialDerivativeTimbreScale(t[:,0],t[:,1],s,i))
  82.     return df
  83.  
  84. #MODIFIED TO ACCEPT ARGS (f, a, s, i)
  85. #Calculates an individual partial derivative of the dissonance curve
  86. cpdef DTYPE_t partialDerivativeTimbreScale(np.ndarray[DTYPE_t, ndim=1] f, np.ndarray[DTYPE_t, ndim=1] a, np.ndarray[DTYPE_t,ndim=1] s, unsigned int differentiatedFreq):
  87.     length = f.size
  88.     assert differentiatedFreq < length
  89.     cdef np.ndarray[DTYPE_t, ndim=1] d
  90.     d = np.array([])
  91.     cdef unsigned int i
  92.     for i in range(length):
  93.         if i == differentiatedFreq:
  94.             d = np.append(d,1)
  95.         else:
  96.             d = np.append(d,0)
  97.  
  98.     # Up to this point, all we need to go is add up all of the grad1Timbre( union([f, a, d], [f*s[i], a, d*s[i]]) )
  99.     cdef DTYPE_t currGrad
  100.     currGrad = 0.0
  101.     for i in range(s.size):
  102.         currGrad += grad1Timbre(np.append(f, f*s[i]), np.append(a, a), np.append(d, d*s[i]) )
  103.     return currGrad
  104.  
  105. #MODIFIED TO ACCEPT 3 1D INPUTS, MODIFIED TO ACCEPT UNSORTED
  106. #A component of partialDerivativeTimbreScale. Not all that useful by itself.
  107. cdef DTYPE_t grad1Timbre(np.ndarray[DTYPE_t,ndim=1] f, np.ndarray[DTYPE_t,ndim=1] a, np.ndarray[DTYPE_t,ndim=1] d):
  108.     cdef DTYPE_t runningGrad1
  109.     runningGrad1 = 0.0
  110.     cdef unsigned int length = f.size
  111.     cdef unsigned int i
  112.     cdef unsigned int j
  113.     cdef unsigned int ii
  114.     cdef unsigned int jj
  115.    
  116.     # v Below: For Optimization Attempt
  117.     cdef DTYPE_t sbase
  118.     cdef DTYPE_t sbasesq
  119.     cdef DTYPE_t df
  120.     # ^ Above: For Optimization attempt
  121.    
  122.     for i in range(length):
  123.         for j in range(i+1,length):
  124.             if not( d[i] < 0.5 and d[j] < 0.5 ):
  125.                
  126.                 if f[j] >= f[i]:
  127.                     (ii,jj) = (i,j)
  128.                 else:
  129.                     (ii,jj) = (j,i)
  130.                
  131.                 if (d[ii] > .5) and (d[jj] < .5):
  132.                     # runningGrad1 += float_min(float_abs(a[ii]),float_abs(a[jj]))*gradF1(f[ii],f[jj])*d[ii]
  133.                    
  134.                     # Below: For Optimization Attempt
  135.                     sbase = 1./(0.021*f[ii] + 19)
  136.                     sbasesq = sbase*sbase
  137.                     df = f[jj]-f[ii]
  138.                     runningGrad1 += float_min(float_abs(a[ii]),float_abs(a[jj]))*(cExp(-0.84*df*sbase)*(0.01764*df*sbasesq + 0.84*sbase) - cExp(-1.38*df*sbase)*(0.02898*df*sbasesq + 1.38*sbase))*d[ii]
  139.                 elif (d[ii] < .5) and (d[jj] > .5):
  140.                     # runningGrad1 += float_min(float_abs(a[ii]),float_abs(a[jj]))*gradF2(f[ii],f[jj])*d[jj]
  141.                    
  142.                     # Below: For Optimization Attempt
  143.                     sbase = 1./(0.021*f[ii] + 19)
  144.                     df = f[jj]-f[ii]
  145.                     runningGrad1 += float_min(float_abs(a[ii]),float_abs(a[jj]))*(1.38*sbase*cExp(-1.38*df*sbase) - 0.84*sbase*cExp(-0.84*df*sbase))*d[jj]
  146.                 else:
  147.                     # runningGrad1 += float_min(float_abs(a[ii]),float_abs(a[jj]))*(d[ii]*gradF1(f[ii],f[jj]) + d[jj]*gradF2(f[ii],f[jj]))
  148.                    
  149.                     # Below: For Optimization Attempt
  150.                     sbase = 1./(0.021*f[ii] + 19)
  151.                     sbasesq = sbase*sbase
  152.                     df = f[jj]-f[ii]
  153.                     runningGrad1 += float_min(float_abs(a[ii]),float_abs(a[jj]))*(d[ii]*(cExp(-0.84*df*sbase)*(0.01764*df*sbasesq + 0.84*sbase) - cExp(-1.38*df*sbase)*(0.02898*df*sbasesq + 1.38*sbase)) + d[jj]*(1.38*sbase*cExp(-1.38*df*sbase) - 0.84*sbase*cExp(-0.84*df*sbase)))
  154.     return runningGrad1
  155.  
  156. # The following code is for calculating the amplitude gradients.
  157.  
  158. # Calculates the gradient of the dissonance curve (does not consider amplitudes)
  159. cpdef np.ndarray[DTYPE_t, ndim=1] aGradTimbreScale(np.ndarray[DTYPE_t, ndim=2] t,np.ndarray[DTYPE_t,ndim=1] s):
  160.     cdef np.ndarray[DTYPE_t, ndim=1] df
  161.     df = np.array([])
  162.     cdef unsigned int i
  163.     for i in range(0, t.shape[0]):
  164.         df = np.append(df,aPartialDerivativeTimbreScale(t[:,0],t[:,1],s,i))
  165.     return df
  166.  
  167. #Calculates the amplitude partial derivatives of the dissonance curve.
  168. cpdef DTYPE_t aPartialDerivativeTimbreScale(np.ndarray[DTYPE_t, ndim=1] f, np.ndarray[DTYPE_t, ndim=1] a, np.ndarray[DTYPE_t,ndim=1] s, unsigned int differentiatedAmp):
  169.     length = f.size
  170.     assert differentiatedAmp < length
  171.     cdef np.ndarray[DTYPE_t, ndim=1] d
  172.     d = np.array([])
  173.     cdef unsigned int i
  174.     for i in range(length):
  175.         if i == differentiatedAmp:
  176.             d = np.append(d,1)
  177.         else:
  178.             d = np.append(d,0)
  179.  
  180.     # Up to this point, all we need to go is add up all of the aGrad1Timbre( union([f, a, d], [f*s[i], a, d) )
  181.     cdef DTYPE_t currGrad
  182.     currGrad = 0.0
  183.     for i in range(s.size):
  184.         currGrad += aGrad1Timbre(np.append(f, f*s[i]), np.append(a, a), np.append(d, d) )
  185.     return currGrad
  186.  
  187. #A component of aPartialDerivativeTimbreScale. Not all that useful by itself.
  188. cdef DTYPE_t aGrad1Timbre(np.ndarray[DTYPE_t,ndim=1] f, np.ndarray[DTYPE_t,ndim=1] a, np.ndarray[DTYPE_t,ndim=1] d):
  189.     cdef DTYPE_t runningGrad1
  190.     runningGrad1 = 0.0
  191.     cdef unsigned int length = f.size
  192.     cdef unsigned int i
  193.     cdef unsigned int j
  194.     cdef unsigned int ii
  195.     cdef unsigned int jj
  196.    
  197.     # v Below: For Optimization Attempt
  198.     cdef DTYPE_t s
  199.     # ^ Above: For Optimization attempt
  200.    
  201.     for i in range(length):
  202.         for j in range(i+1,length):
  203.             if not( d[i] < 0.5 and d[j] < 0.5 ):
  204.                
  205.                 if f[j] >= f[i]:
  206.                     (ii,jj) = (i,j)
  207.                 else:
  208.                     (ii,jj) = (j,i)
  209.                
  210.                 if ( d[ii] < 0.5 and d[jj] > 0.5 ):
  211.                     s = 0.24/(0.021*f[ii] + 19)
  212.                     runningGrad1 += float_dMinAbs(a[jj], a[ii]) * (cExp(-3.5*s*(f[jj]-f[ii])) - cExp(-5.75*s*(f[jj]-f[ii])) )
  213.                
  214.                 elif( d[ii] > 0.5 and d[jj] < 0.5):
  215.                     s = 0.24/(0.021*f[ii] + 19)
  216.                     runningGrad1 += float_dMinAbs(a[ii], a[jj]) * (cExp(-3.5*s*(f[jj]-f[ii])) - cExp(-5.75*s*(f[jj]-f[ii])) )
  217.                
  218.                 else:
  219.                     assert a[ii] == a[jj]
  220.                     s = 0.24/(0.021*f[ii] + 19)
  221.                     runningGrad1 += float_sign(a[ii]) * (cExp(-3.5*s*(f[jj]-f[ii])) - cExp(-5.75*s*(f[jj]-f[ii])) )
  222.     return runningGrad1
  223.  
  224.  
  225.  
  226. cdef inline DTYPE_t float_dMinAbs(DTYPE_t diffd, DTYPE_t notDiffd): return 0 if float_abs(notDiffd) < float_abs(diffd) else float_sign(diffd)
  227.  
  228. cdef DTYPE_t float_sign(DTYPE_t t):
  229.     if t > 0: return 1
  230.     elif t < 0: return -1
  231.     else: return 0
  232.  
  233. # The following is solely for external use:
  234.  
  235. # Calculates the dissonance of two partials of the form [f,a]
  236. cpdef DTYPE_t diss2Partials(np.ndarray[DTYPE_t,ndim=1] p1, np.ndarray[DTYPE_t,ndim=1] p2):
  237.     cdef DTYPE_t f1 = p1[0]
  238.     cdef DTYPE_t f2 = p2[0]
  239.     cdef DTYPE_t a1 = float_abs(p1[1])
  240.     cdef DTYPE_t a2 = float_abs(p2[1])
  241.    
  242.     cdef DTYPE_t s
  243.     s = 0.24/(0.021*f1 + 19)
  244.     return (float_min(a1,a2) * (cExp(-3.5*s*(f2-f1)) - cExp(-5.75*s*(f2-f1)) ) )
  245.  
  246. # NEEDS FIXING
  247. # Calculates the dissonance of two timbres
  248. cpdef DTYPE_t diss2Timbres(np.ndarray[DTYPE_t,ndim=2] t1, np.ndarray[DTYPE_t,ndim=2] t2):
  249.     return diss1Timbre(np.append(t1[:,0],t2[:,0]), np.append(t1[:,1],t2[:,1]) )
  250.  
  251. #Modified to reflect adjusted derivative-multiplier for gradmatrices.
  252. cpdef np.ndarray[DTYPE_t,ndim=2] transpose(np.ndarray[DTYPE_t,ndim=2] t, DTYPE_t ratio):
  253.     if t.shape[1] == 2:
  254.         return np.dot(t, np.array([[ratio,0],[0,1]]))
  255.     else:
  256.         return np.dot(t, np.array([[ratio,0,0],[0,1,0],[0,0,ratio]]))
Advertisement
Add Comment
Please, Sign In to add comment
Advertisement
Public Pastes
Advertisement
We use cookies for various purposes including analytics. By continuing to use Pastebin, you agree to our use of cookies as described in the Cookies Policy.  OK, I Understand
Not a member of Pastebin yet?
Sign Up, it unlocks many cool features!