lasp/lasp/lasp_imptube.py

174 lines
5.0 KiB
Python

#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""!
Author: J.A. de Jong - ASCEE
Description: Two-microphone impedance tube methods
"""
__all__ = ['TwoMicImpedanceTube']
# from lrftubes import Air
from .lasp_measurement import Measurement
from numpy import pi, sqrt, exp
import numpy as np
from scipy.interpolate import UnivariateSpline
# from lrftubes import PrsDuct
from functools import lru_cache
class TwoMicImpedanceTube:
def __init__(self, mnormal: Measurement,
mswitched: Measurement,
s: float,
d1: float,
d2: float,
fl: float = None,
fu: float = None,
periodic_method=False,
# mat= Air(),
D_imptube = 50e-3,
thermoviscous = True,
**kwargs):
"""
Initialize two-microphone impedance tube methods
Args:
mnormal: Measurement in normal configuration
mswitched: Measurement in normal configuration
s: Microphone distance
fl: Lower evaluation frequency
fu: Lower evaluation frequency
kwargs: tuple with extra arguments, of which:
N: Period length of periodic excitation *obligatory*
chan0: Measurement channel index of mic 0
chan1: Measurement channel index of mic 1
"""
self.mnormal = mnormal
self.mswitched = mswitched
self.mat = mat
self.s = s
self.d1 = d1
self.d2 = d2
self.fl = fl
if fl is None:
ksmin = 0.1*pi
kmin = ksmin/s
self.fl = kmin*mat.c0/2/pi
self.fu = fu
if fu is None:
ksmax = 0.8*pi
kmax = ksmax/s
self.fu = kmax*mat.c0/2/pi
self.thermoviscous = thermoviscous
self.D_imptube = D_imptube
self.periodic_method = periodic_method
self.channels = [kwargs.pop('chan0', 0), kwargs.pop('chan1', 1)]
# Compute calibration correction
if periodic_method:
self.N = kwargs.pop('N')
freq, C1 = mnormal.periodicCPS(self.N, channels=self.channels)
freq, C2 = mswitched.periodicCPS(self.N, channels=self.channels)
else:
self.nfft = kwargs.pop('nfft', 16000)
freq, C1 = mnormal.CPS(nfft=self.nfft, channels=self.channels)
freq, C2 = mswitched.CPS(nfft=self.nfft, channels=self.channels)
# Store this, as it is often used for just a single sample.
self.C1 = C1
self.freq = freq
self.il = np.where(self.freq<= self.fl)[0][-1]
self.ul = np.where(self.freq > self.fu)[0][0]
# Calibration correction factor
# self.K = 0*self.freq + 1.0
K = sqrt(C2[:,0,1]*C1[:,0,0]/(C2[:,1,1]*C1[:,1,0]))
# self.K = UnivariateSpline(self.freq, K.real)(self.freq) +\
# 1j*UnivariateSpline(self.freq, K.imag)(self.freq)
self.K = K
def cut_to_limits(self, ar):
return ar[self.il:self.ul]
def getFreq(self):
"""
Array of frequencies, cut to limits of validity
"""
return self.cut_to_limits(self.freq)
@lru_cache
def G_AB(self, meas):
if meas is self.mnormal:
C = self.C1
freq = self.freq
else:
if self.periodic_method:
freq, C = meas.periodicCPS(self.N, self.channels)
else:
freq, C = meas.CPS(nfft=self.nfft, channels=self.channels)
# Microphone transfer function
G_AB = self.K*C[:,1,0]/C[:,0,0]
return self.getFreq(), self.cut_to_limits(G_AB)
def k(self, freq):
"""
Wave number, or thermoviscous wave number
"""
if self.thermoviscous:
D = self.D_imptube
S = pi/4*D**2
d = PrsDuct(0, S=S, rh=D/4, cs='circ')
d.mat = self.mat
omg = 2*pi*freq
G, Z = d.GammaZc(omg)
return G
else:
return 2*pi*freq/self.mat.c0
def R(self, meas):
freq, G_AB = self.G_AB(meas)
s = self.s
k = self.k(freq)
d1 = self.d1
RpA = (G_AB - exp(-1j*k*s))/(exp(1j*k*s)-G_AB)
R = RpA*exp(2*1j*k*(s+d1))
return freq, R
def alpha(self, meas):
"""
Acoustic absorption coefficient
"""
freq, R = self.R(meas)
return freq, 1 - np.abs(R)**2
def z(self, meas):
"""
Acoustic impedance at the position of the sample, in front of the sample
"""
freq, R = self.R(meas)
return freq, self.mat.z0*(1+R)/(1-R)
def zs(self, meas):
"""
Sample impedance jump, assuming a cavity behind the sample with
thickness d2
"""
freq, R = self.R(meas)
z0 = self.mat.z0
k = 2*pi*freq/self.mat.c0
d2 = self.d2
zs = 2*z0*(1-R*exp(2*1j*k*d2))/((R-1)*(exp(2*d2*k*1j)-1))
return freq, zs