tonal_audibility implementation

mosqito/sq_metrics/tonality/tonal_audibility/_criteria.py

@@ -0,0 +1,182 @@
+from collections import Counter
+import _critical_band
+import numpy as np
+import operator
+
+#Function that determines if the number of spectral lines over and below ft is at least 5. LS first iteration criteria
+#Param
+#	Inputs:
+#			listIndex - list of index corresponding to the spectral lines that fulfill the LS_ini+6 criteria 
+#			Li_ft_index - index of ft
+#	Output:
+#			OK - boolean that is set to False if the num_lines over or below ft is less than 5, otherwise is set to True
+def _iteration_criteria_LS(listIndex, Li_ft_index):
+	nu=0 #num of lines below ft
+	no=0 #num of lines over de ft
+	for q in range(0, len(listIndex)):
+		if listIndex[q]<Li_ft_index:
+			nu=nu+1
+		if listIndex[q]>Li_ft_index:
+			no=no+1
+	if no<5 or nu<5:
+		print("Less than 5 lines")
+		OK=False
+	else :
+		OK=True
+	return OK
+
+
+#Function that during the iteration for LS' calculation tests if the new mean value is 
+#equal within a tolerance of ±0,005 dB to that of the previous iteration step. LS' second iteration criteria
+#Param: 
+#	Inputs: list_Ls: list that contains the values of LS in every iteration step
+#	Output: stop - boolean that is set to True if the iteration must be stopped because the criteria is met
+def _iteration_criteria_LS_2(list_LS):
+	for j in range(1,len(list_LS)):
+		if abs(list_LS[j]-list_LS[j-1])<0.005:
+			print("dif < de 0.005")
+			stop=True
+		else: 
+			stop=False
+	return stop
+
+
+#The level frequency lines that can be defined as tones, ft, inside a critical band: a tone may only be present if 
+#the level of the spectral line considered is at least 6 dB greater than the corresponding mean narrow-band level LS.
+#This function tests if the tone satisfies this criteria
+#	Inputs: LS - calculated mean narrow band level for ft
+#			Li - ft's tone level
+#	Output: tone_criteria: boolean set to False if the criteria is not fulfilled 
+def _tone_criteria(LS, Li):
+	if Li>LS+6:
+		tone_criteria=True
+	else: tone_criteria=False
+	return tone_criteria
+
+
+#Function that determines if the tone meets the distinctness criteria
+#Param: 
+#	Inputs: ft - tone frequency
+#			Lu - tone level of the frequency of the first spectral line below the tone ft
+#			Lo - tone level of the frequency of the first spectral line over the tone ft
+#			fu - tone level of the frequency of the first spectral line below the tone ft
+#			fo tone level of the frequency of the first spectral line over the tone ft
+#			LT_max - maximum narrow-band level (Li) of the tone	
+#			Tone_BW - tone bandwidth: sum of the bandwidths of the spectral lines contributing to the tone
+#	Output: ok - boolean that is set to False when the tone doesn´t meet the distinctness criteria, and to True if it does
+def _distinctness_criteria(ft, Lu,Lo,fu,fo, LT_max, Tone_BW):
+	delta_fR=26*(1+0.001*ft)
+	delta_Lu=(ft/2)*((LT_max-Lu)/(ft-fu))
+	delta_Lo=(ft)*((LT_max-Lo)/(fo-ft))
+	#print("delta_Lu and delta_Lo must be higher or equal to 24 dB. They are: (", round(delta_Lu,2),",", round(delta_Lo,2), ")")
+	#print("The tone bandwidth must be lower than the max: ", round(delta_fR,2))
+	if (Tone_BW<delta_fR) and (delta_Lu>=24) and (delta_Lo>=24):
+		ok=True
+	else: ok=False
+	return ok
+
+
+#Function that checks the dictionary for the tones that are distinct but not audible, and erase them from it
+#	Inputs: dict_lines_assigned - This dictionary is used to see what lines are assigned to each tone in Lt, so they are 
+#			not added twice in Ltm. Keys: ft, values: list of Lis assigned to get ft's Lt
+#			list_ft - list of audible tones
+#			dict_aud - dictionary keys ft; values delta_L (audibilities)
+#			list_not_distinct - list of audible tones that do not meet the distinctness conditions
+#	Output: dict_lines_assigned - updated dictionary of which lines are assigned to each tone in Lt, without the discarded
+#			tones (keys) that are not distinct if any
+#			list_ft - updated list of audible tones without the discarded tones that are not distinct if any
+#			dict_aud -updated dictionary of audibilities, without the discarded tones (keys) that are not distinct if any
+def _discard_not_distinct(dict_lines_assigned, list_ft, dict_aud,list_not_distinct):
+	dict_lines_assigned_aux=dict_lines_assigned.copy()
+	#Erase tones, that are finally not present, from the dictionary (distinct but not audible/present)
+	for keys in dict_lines_assigned_aux:
+		if(keys in list_ft)==False:
+			del dict_lines_assigned[keys] 
+
+	#Erase from the list of audible tones those that are not supposed to be studied even though they are audible 
+	#I.e. those which do not meet distinctness criteria
+	for z in range(0, len(list_not_distinct)):
+		f=list_not_distinct[z]
+		if f in list_ft:
+			list_ft.remove(f) 
+		if f in dict_aud:
+			del dict_aud[f]
+	return dict_lines_assigned, list_ft, dict_aud
+
+
+#"It is possible for the energy of individual spectral lines to be assigned to a number of neighbouring tones
+#at the same time. Upon addition of the tone levels of neighbouring tones, the energy of these individual
+#spectral lines may not be summed more than once."
+#This function ensures that when summing the tone levels of neighbouring tones, the energy of these individual spectral lines 
+#are not summed more than once. Regarding the calculation of Ltm.
+#Param: 
+#	Inputs: dict_lines_assigned - This dictionary is used to see what lines are assigned to each tone in Lt, so they are 
+#			not added twice in Ltm. Keys: ft, values: list of Lis assigned to get ft's Lt
+#			dict_aud - dictionary keys ft; values delta_L (audibilities)
+#			list_ft_audib - list of audible tonos in the spectra (yet to be assessed whether they meet all the criteria)
+#			list_freqs - list of all the frequencies in the spectra, neccesary to get the index of the tones that
+#						use the same level(s) to obtain Lt
+#			list_Li - list of all the levels of respective f in the spectra, once obtained the index, the level of the tones
+#					is determined with this list. If one is to be erased it would be the least pronounced tone
+#	Output: dict_lines_assigned - updated dictionary, without the deleted tones (keys) if any
+#			list_ft_audib_aux - updated list of audible tones, without the deleted tones (keys) if any
+#			dict_aud - updated dictionary of audibilities, without the deleted tones (keys) if any
+def _single_addition_lines(dict_lines_assigned, dict_aud, list_ft_audib, list_freqs, list_Li):
+	list_aux1=[]
+	list_aux2=[]
+	list_ft_audib_aux=list_ft_audib.copy()
+	for m in range(1, len(list_ft_audib)-1):
+		key1=list_ft_audib[m-1]
+		key2=list_ft_audib[m]
+		#Get list of Lis used to get each of the tones Lt
+		list_aux1=dict_lines_assigned.get(key1)
+		list_aux2=dict_lines_assigned.get(key2)
+		index_tone1=list_freqs.index(key1)
+		index_tone2=list_freqs.index(key2)
+		Li_tone1=list_Li[index_tone1]
+		Li_tone2=list_Li[index_tone2]
+		#Counter is used to Compare if unordered lists have the same elements
+		if(Counter(list_aux1)==Counter(list_aux2))==True:
+			#In case they have the same lines, keep only the most pronounced tone
+			if(Li_tone1>Li_tone2):
+				list_ft_audib_aux.remove(list_ft_audib[m])
+				del dict_lines_assigned[key2]
+				del dict_aud[key2]
+			else : 
+				list_ft_audib_aux.remove(list_ft_audib[m-1])
+				del dict_lines_assigned[key1]
+				del dict_aud[key1]
+	return dict_lines_assigned, list_ft_audib_aux, dict_aud
+
+
+#"If exactly 2 tones with tone frequencies, fT1 and fT2, appear in one critical band, then they are evaluated
+#separately if both tone frequencies lie below 1 000 Hz and the frequency difference, fD."
+#This function is used to check whether two tones in the same CB should be evaluated separately, according to the standard.
+#Param: 
+#	Inputs: list_ft_audib - list of audible tonos in the spectra (after checking all criteria)
+#			dict_aud - dictionary keys ft; values delta_L (audibilities, after checking all criteria)
+#			list_freqs - list of all the frequencies in the spectra, necessary to determine the CB, f1 and f2,
+#						and which tones lie within a CB
+#	Output: list_ft_audib - updated list of audible tones, if two tones can´t be evaluated separately the lowest of those
+#							is erased from the list
+#			dict_aud_aux - updated dictionary of audibilities
+def _aud_tones_within_critical_band(list_ft_audib, dict_aud, list_freqs):
+	dict_aud_aux=dict_aud.copy()
+	for f in dict_aud:
+		dict_aud_CB={}
+		delta_fc,f1,f2=_critical_band._critical_band(f, list_freqs)
+		#check which of the possible tones to be evaluated are contained in the same critical band
+		for f_key in dict_aud:
+			if f_key>f1 and f_key<f2:
+				dict_aud_CB[f_key]=dict_aud.get(f_key)
+		#Get the tone ft with the highest audibility
+		if len(dict_aud_CB)==2:
+			ft=max(dict_aud_CB, key=dict_aud_CB.get)
+			if ft>50 and ft<1000:
+				fd=21*np.power(10,1.2*np.power(abs(np.log10(ft/212)),1.8))
+				#to be evaluated separately both tones must be under 1000Hz and fd>(ft-f)
+				for f_key2 in dict_aud_CB: #check conditions for the tones within ft's CB
+					if (f_key2<1000) and (ft!=f_key2) and (ft-f_key2<fd) and (f_key2 in dict_aud_aux):
+						del dict_aud_aux[f_key2]
+						list_ft_audib.remove(f_key2)
+	return list_ft_audib, dict_aud_aux

mosqito/sq_metrics/tonality/tonal_audibility/_critical_band.py

@@ -0,0 +1,35 @@
+
+import numpy as np
+
+#Function that determines the width and limiting frequencies of the critical band with centre at ft
+#(considering the line spacing for the values of f1 and f2)
+#ISO 20065 5.2 Width Δfc of the critical band (Formulas 2,4,5)
+#Param:
+#	Inputs: ft - frequency of the tone under study, center of the CB
+#			total_list_f - complete list of frequencies of the spectra to determine which freqs lie within the CB
+#	Output: delta_fc - Bandwidth of the critical band
+#			f1 and f2 - limit frequencies of teh CB
+def _critical_band(ft, total_list_f):
+	delta_fc = 25 + 75 * np.power(1 + 1.4 * np.power(ft/ 1000, 2), 0.69)
+	f1=(-delta_fc/2)+(np.power((np.power(delta_fc,2)+4*np.power(ft,2)),1/2))/2
+	f2=f1+delta_fc
+	if f1<total_list_f[0]:
+		f1_=total_list_f[0]
+	else:
+		for i in range(0,len(total_list_f)):
+			if f1>=total_list_f[i]:
+				f1_=total_list_f[i]
+
+	max=len(total_list_f)-1	
+	if f2>total_list_f[max]:
+		f2_=total_list_f[max]
+	else:
+		for m in range(0, len(total_list_f)):
+			if f2>=total_list_f[m]:
+				f2_=total_list_f[m]
+	return delta_fc,f1_,f2_
+
+
+
+
+

mosqito/sq_metrics/tonality/tonal_audibility/_formula.py

@@ -0,0 +1,159 @@
+import numpy as np
+import math
+
+#Function to delimit the line spacing. It returns the N (number of points for the fft) neccesary to get that delta_f.
+#The value of delta_f will be recalculated again afterwards when the number of points for the fft is rounded 
+#to a power of 2
+def _delta_f(fs):
+	#>1.9Hz <4Hz
+	delta_f=3 #duration(of each clip)=1/delta_f
+	N=fs/delta_f;
+	return N,delta_f
+
+
+#Function that returns the smallest power of two that is greater than or equal to the absolute value of N.
+#It also recalculates delta_f according to that new N (number of points for the fft)
+#Param: N - numbre 
+def _nextpow2(N,fs):
+	n=np.log2(N)
+	if(n%2!=0):
+		n=math.ceil(n)
+	else: n=n
+	print(n)
+	new_N=np.power(2,n)
+	#recalculate delta_f
+	delta_f=fs/new_N
+	print(delta_f)
+	return new_N,delta_f
+
+
+#Function that calculates the masking index
+#ISO 20065 5.3.6. Masking index (Formula 13)
+#Param: 
+#	Inputs: ft - tone frequency: frequency of the spectral line (or mid-band frequency of the narrow-band filter), to the level of
+#								 which the tone contributes most strongly
+#	Output: av - masking index: audibility threshold for a specific sound in the presence of a masking sound
+def _masking_index(ft):
+	av=-2-np.log10(1+(np.power(ft/502, 2.5)))
+	return av
+
+
+#Function that calculates the critical band level
+#ISO 20065 5.3.5 Critical band level, LG, of the masking noise (Formula 12)
+#Param: 
+#	Inputs: LS - mean narrow band level: energy mean value of all narrow-band levels in a critical band that does not exceed
+#										 this mean value by more than 6 dB
+#			delta_f - line spacing: distance between neighbouring spectral lines
+#			delta_fc - bandwith of the critical band
+#	Output: LG - critical band level: level of noise that is assigned to the critical band that describes the masking 
+#									  characteristic of the noise for one or more tones of the noise in this critical band
+def _critical_band_level(LS, delta_f, delta_fc):
+	LG=LS+(10*np.log10(delta_fc/delta_f))
+	return LG
+
+
+#Function that calculates the audibility
+#ISO 20065 5.3.7 Determination of the audibility, ΔL (Formula 14)
+#Param: 
+#	Inputs: LT - tone level: energy summation of the narrow-band level with the tone frequency, fT, and the lateral lines
+#							 about fT, assignable to this tone
+#			LG - critical band level: level of noise that is assigned to the critical band that describes the masking 
+#			av - masking index: audibility threshold for a specific sound in the presence of a masking sound
+#	Output: delta_L - audibility: difference between the tone level and the masking threshold				
+def _audibility(LT, LG, av):
+	delta_L=(LT-LG-av)
+	return delta_L
+
+
+
+#Function that calculates LS 
+#ISO 20065 5.3.2 Determination of the mean narrow-band level LS of the masking noise (Formula 6)
+#Param: 
+#	Inputs: M - number of spectral lines except the one under study:
+#			listTones_Li - list of the levels Li of the tones within the CB except for ft
+#			delta_f - line spacing: distance between neighbouring spectral lines
+#			delta_fe - effective bandwidth: is the effective bandwidth in Hz; if a Hanning window is used then the effective bandwidth, Δfe,
+#											is 1,5 times the frequency resolution (line spacing), Δf.
+#	Output: LS - mean narrow band level after the first iteration. Initial LS
+def _mean_narrow_band_level(M, listTones_ini, delta_f, delta_fe):
+	Sum=0
+	for i in range(0,len(listTones_ini)):
+		Sum=Sum+(np.power(10,0.1*listTones_ini[i]))
+	LS=(10*np.log10((1/M)*Sum))+(10*np.log10(1/delta_fe))
+	return LS
+
+
+#Function that calculates the tone level
+#ISO 20065 5.3.3 Determination of the tone level LT of a tone in a critical band (Formula 8)
+#Param: 
+#	Inputs: listTones - list of the levels Li of the tones within the CB that are under LS+6dB
+#			delta_f - line spacing: distance between neighbouring spectral lines
+#			LS - mean narrow band level: energy mean value of all narrow-band levels in a critical band that does not exceed
+#										 this mean value by more than 6 dB
+#			ft_index - index of the tone under study
+#			delta_fe - effective bandwidth: is the effective bandwidth in Hz; if a Hanning window is used then the effective bandwidth, Δfe,
+#											is 1,5 times the frequency resolution (line spacing), Δf.
+#	Output: LT - tone level: energy summation of the narrow-band level with the tone frequency, fT, and the lateral lines
+#							 about fT, assignable to this tone
+#			Tone_BW - tone bandwidth: sum of the bandwidths of the spectral lines contributing to the tone
+#			LT_max - maximum narrow-band level (Li) of the tone			
+def _tone_level(listTones, delta_f, LS, ft, Li, ft_index, delta_fe, dict_lines_assigned):
+	list_lines_assigned=[]
+	list_lines_assigned.append(Li)
+	Sum=np.power(10,0.1*Li) #Add ft's Li to the summatory 
+	LT_max=Li #for distinctness criteria
+	Tone_BW=delta_f #Tone bandwith necessasry for distinctness criteria
+	index_aux=ft_index
+	#below the tone
+	for j in range(ft_index-1,0,-1):
+		if (listTones[j]>Li-10) and (listTones[j]>LS+6):
+			if (j==index_aux-1): 
+				index_aux=j
+				Sum=Sum+(np.power(10,0.1*listTones[j]))
+				list_lines_assigned.append(listTones[j])
+				Tone_BW=Tone_BW+delta_f
+	index_aux=ft_index
+	#over the tone
+	for x in range(ft_index+1,len(listTones)):
+		if (listTones[x]>Li-10) and (listTones[x]>LS+6): 
+			if (x==index_aux+1):
+				index_aux=x
+				list_lines_assigned.append(listTones[x])
+				Sum=Sum+(np.power(10,0.1*listTones[x]))
+				Tone_BW=Tone_BW+delta_f
+	dict_lines_assigned[ft]=list_lines_assigned
+	if Tone_BW>delta_f: 
+		LT=(10*np.log10(Sum))+(10*np.log10(1/delta_fe))
+	else: #if only ft's Li is in the summation
+		LT=10*np.log10(Sum)
+	return LT, Tone_BW, LT_max, dict_lines_assigned
+
+
+#Function that adds the tone level (LT) of all the audible tones within ftm's critical band
+#ISO 20065 5.3.8. Determination of the decisive audibility, ΔLj, of a narrow-band spectrum -Step 3 (Formula 17)
+def _sum_tone_level(list_LT, list_Li):
+	Sum=0
+	for z in range(0,len(list_Li)):
+		Sum=Sum+(np.power(10,0.1*list_Li[z]))
+	if len(list_Li)>2*len(list_LT):
+		LTm=10*np.log10(Sum)-1.76
+	else:
+		LTm=10*np.log10(Sum)
+	return LTm
+
+
+#Function for the determination of the mean audibility of a number of spectra
+#ISO 20065 5.3.9 Determination of the mean audibility ΔL of a number of spectra (Formula 20)
+#The decisive audibility ΔLj is calculated for each narrow-band averaged spectrum (run
+#index j, J is the number). These J audibilities are averaged in energy terms to yield the mean audibility
+#Param:
+#	Inputs: list_decisve_aud - decisive audibility of each spectrum, to be averaged
+#	Output: mean_aud - mean audibility
+def _mean_aud(list_decisive_aud):
+	Sum=0
+	J=len(list_decisive_aud)
+	for i in range(0,len(list_decisive_aud)):
+		Sum=Sum+np.power(10,0.1*list_decisive_aud[i])
+	mean_aud=10*np.log10((1/J)*Sum)
+	return mean_aud
+