Calculs audio numérique pour les nuls

[alka]

C'est vrai qu'on utilise des types de filtres connus dans DStudio. Avec ce qui a été fait, on a déjà quelques valeurs. Mais des tables quelle que soit Fc, c'est long.

Le passage en temporel est peut etre une option simple. On mesure le délai entre les deux pics d'impulsion. Avec le nouveau DStudio v5 on aura les outils. A suivre

[alka]

Encart sur le RTA au 1/n d'octave, bruit blanc et rose

Explication un peu approximative et que j'espère claire :

Le spectre affiche le SPL par fréquence. Pour un bruit blanc, la courbe SPL est plate. Pour un bruit rose, c'est une pente descendante de 3dB/o

Bruit blanc et rose issus de REW :


Le "RTA au 1/n octave" montre une chose un peu différence : affiche l'energie contenue dans des bandes de fréquence découpées par nième d'octave.

Pourquoi faire ces analyses de spectre en "RTA au 1/n d'octave" ? parce que ce type d'analyse correspond plus justement à la perception de l'oreille humaine.

Le découpage en nième d'octave coupe chaque octave en un même nombre n bandes. n est parfois appelé octaveband.
Les bornes de ces bandes sont Fi tel que : Fi = Fi-1 * 2 ^(1/octaveband) ]. F0 est la fréquence de départ.
Ca peut aussi s'écrire Fi = F0 * 2^(i/OctaveBand).

La largeur des bandes augmente proportionnellement à la fréquence Fi.
Le RTA est toujours affiché avec une échelle des fréquences logarithmique ce qui a pour effet de montrer graphiquement des barres de largeur fixe. La largeur des barres dépend de octaveband.

Sur ce type d'analyse, un bruit blanc est montant de 3dB/o et un bruit rose est plat.

Le même bruit rose que précédemment :


Calcul:
Pour la freq bande [Fi, Fi-1] on affiche l'energie du signal contenue dans cette bande.
On peut imager "l'energie contenue dans une bande de freq" comme étant la multiplication (niveau moyen du signal dans la bande) * (largeur de la bande).

Un bruit blanc possede le même niveau a toutes les fréquences. L'energie contenue dans les bandes est croissante. En effet, (niveau constant * largeur de bande qui grandit) augmente. Le RTA d'un bruit blanc sera montant de 3dB/octave.

Un bruit rose a, par construction, un niveau qui baisse régulièrement avec les fréq (baisse de 3dB/octave)
Le RTA de ce signal est plat. En effet, la baisse de niveau du signal compense pile l'augmentation de largeur des bandes de fréquences nème octave.

La formule précise pour l'affichage en RTA est :
RTA (Fi) = 20 * log[ (moyenne modules entre Fi et Fi-1) * racine (Fi - Fi-1)]
tiens on dirait la mesure habituelle de bruit en Vrms par racine de hertz, convertie en dB
Dans la formule, le module est le module de la valeur complexe issu de la FFT.


Le calcul de la "moyenne des modules entre Fi-1 et Fi" se fait soit par moyenne géométrique (on fait la somme et on divise par le nombre) soit en rms (racine de la somme des carrés divisée par le nombre). rms est le plus habituel et plus correct si signal très variable.

exemple :
a gauche issu de la table FFT, les freq en echelle linéaire et le module correspondant
a droite, un extrait de la table de fréq au 1/3 d'octave en partant de 10Hz

4 valeurs calculées pour l'exemple
détail de celle pour 20Hz surlignée en jaune :
la bande va de Fi-1 = 15,9 à Fi=20Hz. largeur = 20 - 15,9 = 4,1
valeur RMS des modules dans cette bande = sqrt((56,4²+54,1²+52,1²)/3) = 54,23
RTA = 20 * log (54,23 * sqrt(4,1) ) = 40,8

[thierryvalk]

Exactement ce qu'il me fallait et en plus super bien expliqué, y a plus qu'a coder cela dans DStudio.
Un tout grand merci

[jimbee]

C'est possible de calculer le délai a appliquer depuis les coefficients des biquads du filtre passe-bas ?

Bonjour,

Petit retour, sauf erreur, pour chaque ordre 2 défini par Q et F , le délai à appliquer est donné par : 1/ ( Q . 2 . Pi . F ) soit 1/ (Q . oméga ) ,
et pour l'ordre 1, c'est 1/ oméga.

crd

[alka]

oula c'est du déterrage.
j'avais oublié le contexte de la question... c'était pour le filtrage soustractif.
ici viewtopic.php?p=4219#p4219 et suivants

pour l'exemple du LR24 @ 1kHz:
LR24 = 2 * but ordre 2 en cascade , chaque Q=0,707

Sur la page le ton blog, http://jimbee.over-blog.com/page-2554286.html
le délai théorique annoncé était de 45ms. Ma simulation de l'époque donnait 0,50ms

La formule pour un ordre 2 : 1/ ( Q . 2 . Pi . F ) = 1/(0,707*6,28*1000) = 0,224ms
deux fois c'est 0,45ms. cqfd.

C'est super d'avoir une formule aussi simple de calcul direct du délai.
merci d'etre revenu pour le signaler
j'ai vu que tu as mis ton blog a jour aussi.

[thierryvalk]

J'ai lu ici :
https://www.pitt-pladdy.com/blog/_20120 ... st_Tools_/
Que sur un sweep, il était facile d'extraire les harmoniques.

C'est encore du domaine des maths.

J'ai de morceaux de code, dont lsconv qui nous sert à produire l'impulse. Si c'est vraiment pas si compliqué, ce serait sympa d'avoir H2 et H3 sur ce signal.

[alka]

j'ai regardé en diagonale.
conaissant précisément le sweep, on doit pouvoir le trouver dans la réponse, l'extraire et regarder les harmoniques mais bon...
"The analyser runs through the file it is given and runs a windowed FFT on it section by section to identify the sweep. It uses various algorithms to work out the rate (assuming an exponential sweep to match normal log frequency scales) and position of the sweep in the file while detecting and ignoring noisy or potentially unreliable sections."

mwouais... plus facile a dire qu'a faire

[thierryvalk]

J'ai commencé à regarder le code source de DRC qui fait ce boulot.
On y parle aussi d'harmonique, mais pas simple car si le code source est documenté, c'est en italien.

Je pense que la disto (importante) est un élément à prendre en compte lors des corrections car s'il y a disto, il y a modification du SPL plus loin.
Et corriger un SPL issu d'une disto ne ferais que du mal.

D'une part on ne corrigerait rien, mais en plus en réduisant le SPL cette disto serait moins noyée dans le signal.

[thierry38]

Hello Thierry,
Il faudrait utilisé la méthode Farina (HOLM,REW).

elle permet d'extraire l'impulse linéaire (càd composée uniquement de pures sinusoides).
c'est l'impulse principale,dégagée de ses harmoniques.(harmoniques que l'on retrouve en amont (temps) sur l'impulse ).

Mais bon,même 1% de disto sur un sinus,la valeur efficace ne varie pas significativement.
Par contre le méthode Farina est séduisante,techniquement parlant (et surtout très rapide).

http://www.nssc.re.kr/?module=file&act=procFileDownload&file_srl=4083&sid=2fb23dbcce0ca1db2cd4e2999f3034b5
il y a plein de doc technique sur le site d'AudioPrecision

[thierryvalk]

Merci pour le lien.
J'avais vu la méthode Farina, mais là c'est très bien expliqué de manière compréhensible par tous.

Sinon, sans trop de rapport, dans mes recherches, je suis tombé sur Acmus:
http://gsd.ime.usp.br/acmus/english/softe.html

Pas encore eu le temps de jouer avec, mais regroupe un ensemble d'outils pour l'analyse du son dans une salle.
Basé sur Eclipse, Opensource en Java.

[thierryvalk]

Il est toujours intéressant de comprendre ce que l’on utilise.
Le Sweep pour les mesures suivant la méthode du Professeur Farina. Voici ce que j’en ai compris, merci de me corriger et compléter.
Tous les détails :
http://pcfarina.eng.unipr.it/Public/Pap ... AES122.pdf

On commence par créer un fichier avec le sweep exponentiel. Les paramètres usuels sont la durée du Sweep, fréquence de départ et fréquence de fin.
Lors de cette création, sera créé un second fichier qui est son inverse.
Ensuite le fichier Sweep sera joué et enregistré.

Pour l’enregistrement, 2 options :
Soit on n’enregistre que le signal de mesure ou l’on enregistre également un signal de référence.
Sais pas trop ce que REW fait avec, c’est pour une référence de timing alors que d’autres système tels que Lsconv que l’on utilise avec DStudio vont utiliser ce signal comme correcteurs d’amplitudes et de phase.

Déjà là, on voit que si l’on veut mesurer un HP de manière précise, on ne prendra pas ce signal sur l’autre canal de la carte son, mais repiqué en sortie d’ampli.
Sortie de l’ampli, un réseau de résistance, une protection éventuelle et l’on rentre dans la carte son.
Avec cette manip, on va corriger la chaine de mesure au plus proche du HP.

Après la mesure, l’on dispose de 2 ou 3 fichiers.
L’enregistrement.
Le sweep inverse.
Eventuellement un enregistrement de référence.

On effectue une convolution et l’on a notre impulse.
En fait pas tout à fait, on va avoir plusieurs impulses décalées dans le temps.
La dernière, que l’on utilise couramment sera la fondamentale tandis que les précédentes les harmoniques.

Voici avec REW, un sweep sur 256K en 48KHz.
On regarde rarement ce qui se passe avant le T0.

A droite notre IR, et entourées H2 et H3

[thierryvalk]

L’explication est fort simple, lorsque l’on a généré du 100Hz par exemple, on a enregistré du 100Hz mais en même temps du 200Hz, 300Hz…
Vu que l’on repasse en temporel, ces fréquences harmoniques sont bien arrivées avant la fondamentale.

Cela se calcule :
delta t=T.(ln(N)/(ln(Fend/Fstart)
N= ordre de l’harmonique
T= durée du sweep
Fend= fréquence finale
Fstart= fréquence initiale

C’est théorique, en réalité on aurra des valeurs un peu différentes.

En effet, cette méthode surpasse sur beaucoup de points les autres méthodes, mais nécessite un début et une fin.
Et partir d’un silence puis généer du signal plein pot va créer des oscillations sur l’impulse.
Pour limiter l’effet, on va appliquer un fade-in au début du sweep et un fade-out à la fin.
Sorte de fenêtrage sur l’ensemble du sweep.

Donc le tout début n’est pas vraiment fiable, ni l’extrème fin.
On rajoute aussi un petit silence avant de démarrer et à la fin ce qui fait que la durée du sweep réel n’est pas forcément juste.

De ce que j’ai lu, REW applique un Sweep linéaire pour les fréquences inférieures à 10Hz.
Cela dans le but réduire la durée du sweep.

[thierryvalk]

Revenons à nos IR.
Ce qui est intéressants, c’est que l’IR que l’on va utiliser est censée être sans les harmoniques.
Et la réponse en fréquence que l’on va en déduire via une FFT aussi.

On parle de fenêtré ou non-fenêtré, en réalité on fenêtre toujours.
Tout est question de taille de fenêtre et dans certains cas du type de fenêtre.

Par défaut, REW applique un fenêtrage de -125mS à +500mS.
Voici en bleu sur nos impulses, on prend que la fondamentale.

Mais si l'on recule le début de cette fenêtre à -800mS, on va prendre aussi les harmoniques

En rouge le signal affiché avec les valeurs par défaut, en vert avec fenêtre entendue.

[thierryvalk]

Le fenêtrage permettant de prendre une partie du signal, on va déplacer signal et fenêtre pour ne prendre que l’IR précédente.

Ce qui donne en réponse en fréquences :

On aurait du avoir H2, visiblement c'est plus H3, décidément ...

Pour confirmation, notre signal de base en réponse de fréquence avec fenêtrage par défaut et affichage des harmoniques : (oui, j'y ai été un peu fort avec la saturation. )

En comparant les 2 graphes, on voit la fréquence décalée. Normal une H2 à 100Hz, par exemple, donne du 200Hz. Et H3 à 100Hz donne 300Hz.

[thierryvalk]

Mais où est donc cette H2 ?

En toute logique entre H0, la fondamentale, et H3. On zoom :

Pas flagrant.

Petit test en RTA :
HS, REW dispose d'une fonction sur le générateur "Lock frequency to FFT".
Très intelligent, en RTA aussi on ne prend qu'une portion du signal qu'il faut fenêtrer pour limiter les transitions au début et à la fin de cet échantillon.
En collant le signal du générateur à celui de la FFT du RTA, on ne doit plus faire du fenêtrage compliqué, du simple rectangulaire suffit et donne un affichage bien plus fidèle.

Là, on a la confirmation que la H2 (et les autres H paires) est relativement faible.
Pourtant sur le graphe Distorsion de la mesure elle était bien présente.

La cause ? Je présume que mon sweep était trop court, avec une fréquence min trop basse et au final H3 à débordé au niveau temporel sur H2.
Pour vérifier, voici une nouvelle mesure en rallongeant le temps du Sweep tout en commençant à 20Hz et non plus 10Hz.

Pas mal de différences avec le graphe du post précédent.

Il faut donc bien faire attention, la mesure de disto, c'est du signal pas bien élevé et il faut lui laisser un certain "espace" temporel.

Farina, suggère une moyenne de 4 sweeps moyennés pour améliorer les rapport signal bruit sur des durées de 15 à 45 secondes.

[thierryvalk]

Et voici la théorie appliquée à DStudio.

Il y a encore un peu de boulot à déterminer les longueurs des fenêtres, mais l'idée est là.
La prétention n'est pas de faire un système de mesure mais plus de voir s'il se passe quelque chose de pas normal.

[thierryvalk]

A droite de notre du T0 de l'impulse, on a la décroissance du SPL selon le temps.

Sox permet de créer très simplement une image du spectrogramme:

En X le temps, en Y la fréquence et la couleur donne le SPL.
On peut facilement améliorer l'échelle de temps, mais l'échelle de fréquence est linéaire. Ce qui n'est pas très pratique.

Autre représentation, le Waterfall.
Cette fois, en X la fréquence en log comme on a l'habitude, en Y le SPL et en Z il faudrait le temps.
A partir du T0 et en allant vers la droite, on prend des échantillons, fenêtrage puis FFT pour récupérer SPL en fonction de f à un temps donné.

Voici un début dans DStudio :

L'idée est là, mais il faudrait un affichage en 3D avec effet de perspective pour avoir le temps.
Malheureusement la lib pour créer ces graphiques dans DStudio ne le permet pas.
Je me dis qu'il devrait être possible d'avoir une représentation un peu plus correcte.
A mon niveau, je coince surtout sur le fait que l'échelle des fréquences est en log. Pour le niveau, c'est relativement simple, pour chaque courbe on applique un gain au SPL pour dilater l'affichage.
Après, il faudra tracer quelques lignes pour les nouveaux axes, c'est pas trop compliqué.

[thierryvalk]

En dessinant quelques lignes c'est déjà mieux.
Mesure des HP de mon PC.