Source Code for Module wiat.util

  1  # coding: latin-1 
  2  # 
  3  # Copyright 2008 Antoine Lefebvre 
  4  # 
  5  # This file is part of "The Wind Instrument Acoustic Toolkit". 
  6  # 
  7  # "The Wind Instrument Acoustic Toolkit" is free software: 
  8  # you can redistribute it and/or modify it under the terms of 
  9  # the GNU General Public License as published by the 
 10  # Free Software Foundation, either version 3 of the License, or 
 11  # (at your option) any later version. 
 12  # 
 13  # "The Wind Instrument Acoustic Toolkit" is distributed in the hope 
 14  # that it will be useful, but WITHOUT ANY WARRANTY; without even 
 15  # the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. 
 16  # See the GNU General Public License for more details. 
 17  # 
 18  # You should have received a copy of the GNU General Public License 
 19  # along with "The Wind Instrument Acoustic Toolkit". 
 20  # If not, see <http://www.gnu.org/licenses/>. 
 21  # 
 22  # All rights reserved. 
 23   
 24  """ util.py - classes and functions that are not specific to the calculation 
 25                methods. 
 26   
 27  """ 
 28   
 29  import numpy 
 30  import scipy 
 31   
 32  import wiat 
 33   
 34 -def air_properties(T): 
 35      """ Properties of air as a function of temperature 
 36   
 37          The formula, valid for ±10 Celcius around 26.85 C, comes from Table 1  
 38          of reference 1. 
 39           
 40          References: 
 41           
 42              1. Keefe, D. H. (1984).  
 43                 Acoustical wave propagation in cylindrical ducts: Transmission  
 44                 line parameter approximations for isothermal and nonisothermal  
 45                 boundary conditions. J. Aoucst. Soc. Am., 75(1), 58-62. 
 46       
 47          TODO: look for a more recent source to double check the validity 
 48   
 49      """ 
 50      deltaT = T - 26.85 
 51      properties = { 
 52          # density   (kg/m^3) 
 53          'rho': 1.1769    * (1 - 0.00335*deltaT), 
 54          # viscosity (kg/m sec) 
 55          'mu': 1.8460E-5 * (1 + 0.00250*deltaT), 
 56          # specific heat ratio 
 57          'gamma': 1.4017    * (1 - 0.00002*deltaT), 
 58          # square root of Prandtl number : sqrt(Cp * mu / kappa) 
 59          'nu': 0.8410    * (1 - 0.00020*deltaT), 
 60          # free space sound speed (m/s) 
 61          'c': 3.4723E+2 * (1 + 0.00166*deltaT)} 
 62      return properties 
 63   
 64 -def calculate_interval_in_cents(f1, f2): 
 65      " Calculate the interval in cents between two frequencies. " 
 66      return 1200. * numpy.log2(f2/f1) 
 67   
 68 -def calculate_ratio_in_db(A1,A2): 
 69      " Calculate the ratio in dB between two values. " 
 70      return 20. * numpy.log10(A2/A1) 
 71   
 72 -def transform_R2Z(R): 
 73      " Convert reflection coefficient to normalized impedance " 
 74      return (1+R)/(1-R) 
 75   
 76 -def transform_Z2R(Z): 
 77      " Convert normalized impedance to reflection coefficient " 
 78      return (Z-1)/(Z+1) 
 79   
 80 -def equivalent_length(R, k): 
 81      """ Return the equivalent length of a complex reflection coefficient  
 82          R: reflection coefficient 
 83          k: wavenumber 
 84      """ 
 85   
 86      #return numpy.real( numpy.log(-R/numpy.abs(R)) / (-1j*2*k) ) 
 87      return numpy.angle(-R)/(-2*k) 
 88   
 89   
 90 -def load_or_calculate(filename, func, **args): 
 91      """ Function that will load a backup of the results if available, 
 92      otherwise it calculates the function and save the results. 
 93      """ 
 94      try: 
 95          f,Z = wiat.Data.load_from_mat(('f', 'Z'), filename) 
 96      except IOError,arg: 
 97          print arg 
 98          print 'No backup available (%s), calculating...'%filename 
 99          Z = func(**args) 
100          scipy.io.savemat(filename, {'f': args['freqs'], 'Z': Z}) 
101      return Z 
102   
103 -def __create_equal_tempered_scale(): 
104      scale = {} 
105      note_name = \ 
106          ['Cx','Dx flat','Dx','Ex flat','Ex',\ 
107           'Fx','Gx flat','Gx','Ax flat','Ax','Bx flat','Bx'] 
108   
109      # create the frequency dictionary for 8 octave 
110      for i in range(8): 
111          for j,note in enumerate(note_name): 
112              name = note.replace("x","%d"%i) 
113              scale[name] = 440*2**(-4+i+(j-9)/12.) 
114      return scale 
115   
116  EQUAL_TEMPERED_SCALE = __create_equal_tempered_scale() 
117